NuSystFitter Class Reference

#include <NuSystFitter.h>

Inheritance diagram for NuSystFitter:

NuABFitter List of all members.

Public Member Functions

 NuSystFitter ()
virtual ~NuSystFitter ()
virtual double Up () const
virtual void SetErrorDef (const Double_t up)
virtual double operator() (const std::vector< double > &pars) const
virtual double Likelihood (const NuMMParameters &mmPars) const
virtual void Setup (const std::vector< std::string > helperInputs, const std::vector< NuHistInterpolator * > fdDataInput, const std::vector< NuHistInterpolator * > ndDataInput, const std::string xsecFile, const std::vector< std::string >, const std::string, Int_t)
virtual void push_back (NuMMRun *run)
virtual const NuMMParameters Minimise (const NuMMParameters &mmPars)
virtual const NuMMParameters Minimise (const NuMMParameters &mmPars, bool &success)
virtual const NuMMParameters FCMinimise (const NuMMParameters &mmPars)
 More robustly minimizes by doing a wide grid search first - same as the FC fitting.
virtual const NuMMParameters FCMinimiseNu (const NuMMParameters &mmPars)
virtual const NuMMParameters FCMinimiseNSI (const NuMMParameters &mmPars)
virtual const TGraph * PlotContours (const NuMMParameters &mmPars)
virtual void NoOscPrediction ()
virtual void NumberOfContourPoints (const Int_t numContourPoints)
virtual const std::pair< Double_t,
Double_t > 
SingleParDm2Errors (const NuMMParameters &mmPars, const bool quiet=false)
virtual const std::pair< Double_t,
Double_t > 
SingleParDm2BarErrors (const NuMMParameters &mmPars, const bool quiet=false)
virtual const std::pair< Double_t,
Double_t > 
SingleParSn2BarErrors (const NuMMParameters &mmPars, const bool quiet=false)
virtual const std::pair< Double_t,
Double_t > 
SingleParErrors (const NuMMParameters &mmPars, const int idx, const bool quiet)
virtual const std::vector<
std::pair< Double_t, Double_t > > 
NeutrinoContourFinder (const NuMMParameters &mmPars)
virtual const std::pair< Double_t,
Double_t > 
DmScanForContour (const NuMMParameters &mmPars, Double_t targetChi2)
virtual void Sn2PointsForNeutrinoContour (const std::vector< Double_t > &sn2PointsForNeutrinoContour)
virtual const TGraph * OneSidedNeutrinoContourFinder (const NuMMParameters &inputPars)
virtual const TGraph * TwoSidedNeutrinoContourFinder (const NuMMParameters &inputPars)
virtual const TGraph * Chi2SliceSin2bar (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minval=0.8, const Double_t maxval=1.2)
virtual const TGraph * Chi2SliceDm2bar (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minval=1e-3, const Double_t maxval=3e-3)
virtual const TGraph * Chi2SliceDm2 (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minval=1e-3, const Double_t maxval=3e-3)
virtual const TGraph * Chi2ValleyBar (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minSin2bar=0., const Double_t maxSin2bar=2.)
virtual const TGraph * Chi2ValleyBarDm2Bar (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minDm2bar=0., const Double_t maxDm2bar=10e-3)
virtual const TGraph * Chi2ValleyDm2 (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minDm2=0., const Double_t maxDm2=10e-3)
virtual const TGraph * Chi2ValleySn2 (const NuMMParameters &mmPars, const Int_t numpoints=100, const Double_t minSn2=0., const Double_t maxSn2=1.)
virtual TH2D * LikelihoodSurfaceBar (const NuMMParameters &mmPars, Int_t sinPoints, Int_t massPoints, Double_t Sn2BarMin=0.0, Double_t Sn2BarMax=1.0, Double_t Dm2BarMin=0.0, Double_t Dm2BarMax=12e-3, Bool_t True_Minimum=true)
 Function to generate a likelihood surface.
virtual Double_t FillLikelihoodSurface (TH2D &surface, NuMMParameters mmPars)
 Fills a histogram with likelihood values. This function accepts an arbitrarily binned histogram, and fills the values according to the likelihood. It returns the value of the minimum likelihood. (This does the CC neutirno spectrum).
virtual Double_t FillLikelihoodSurfaceBar (TH2D &surface, NuMMParameters mmPars)
 Fills a histogram with likelihood values. This function accepts an arbitrarily binned histogram, and fills the values according to the likelihood. It returns the value of the minimum likelihood.
virtual Double_t FillLikelihoodSurfaceNSI (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair)
virtual Double_t FillLikelihoodSurfaceNSIBar (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair)
virtual Double_t FillLikelihoodSurfaceNSIDm2Eps (TH2D &surface, NuMMParameters mmPars, TH2D &sinvals)
virtual Double_t FillLikelihoodSurfaceNSISn2Eps (TH2D &surface, NuMMParameters mmPars, TH2D &massvals)
virtual Double_t FillLikelihoodSurfaceNSIDm2Sn2 (TH2D &surface, NuMMParameters mmPars, TH2D &epsvals)
virtual Double_t FillLikelihoodSurfaceNSI3D (TH3D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceDecay (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair=1)
virtual Double_t FillLikelihoodSurfaceDecoherence (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair=1)
virtual Double_t FillLikelihoodSurfaceDecoherenceBar (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair=1)
virtual Double_t FillLikelihoodSurfaceCPT (TH2D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceSterile (TH2D &surface, NuMMParameters mmPars, Int_t ContourPair, int Ix=1, int Iy=1, int Nx=0, int Ny=0, TH2D *mmPars_surfaces[]=0, TH2D **converge_surface=0, int par_index=-1, int n_iter=1, double *iter=0)
virtual Double_t FillLikelihoodSurfaceBarGMu (TH2D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceGMuGTau (TH2D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceGMuCMuMu (TH2D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceGMuGTau3D (TH3D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceDm2GMuCMu3D (TH3D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransDm2BarGMu (TH2D &surface, NuMMParameters mmPars, TH2D &sinvals)
virtual Double_t FillLikelihoodSurfaceTransDm2GMu (TH2D &surface, NuMMParameters mmPars, TH2D &sinvals)
virtual Double_t FillLikelihoodSurfaceTransDm2GMuBoth (TH2D &surface, NuMMParameters mmPars, TH2D &sinvals)
virtual Double_t FillLikelihoodSurfaceTransSn2BarGMu (TH2D &surface, NuMMParameters mmPars, TH2D &massvals)
virtual Double_t FillLikelihoodSurfaceTransSn2GMu (TH2D &surface, NuMMParameters mmPars, TH2D &massvals)
virtual Double_t FillLikelihoodSurfaceTransSn2GMuBoth (TH2D &surface, NuMMParameters mmPars, TH2D &massvals)
virtual Double_t FillLikelihoodSurfaceTransDm2BarSn2Bar (TH2D &surface, NuMMParameters mmPars, TH2D &gmuvals)
virtual Double_t FillLikelihoodSurfaceTransDm2Sn2 (TH2D &surface, NuMMParameters mmPars, TH2D &gmuvals)
virtual Double_t FillLikelihoodSurfaceTransDm2Sn2Both (TH2D &surface, NuMMParameters mmPars, TH2D &gmuvals)
virtual Double_t FillLikelihoodSurfaceTransGMuBoth (TH1D &surface, TH1D &cmuvals, TH1D &dm2vals, TH1D &sn2vals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransDm2Both (TH1D &surface, TH1D &sn2vals, TH1D &gmuvals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransSn2Both (TH1D &surface, TH1D &dm2vals, TH1D &gmuvals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransCMuMuBoth (TH1D &surface, TH1D &gmuvals, TH1D &dm2vals, TH1D &sn2vals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransCTauTauBoth (TH1D &surface, TH1D &gmuvals, TH1D &dm2vals, TH1D &sn2vals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransGTauBoth (TH1D &surface, TH1D &cmuvals, TH1D &dm2vals, TH1D &sn2vals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransGMuPQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransDm2PQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransSn2PQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransGMuNQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransDm2NQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceTransSn2NQ (TH1D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceDm2BarSn2BarGMu3D (TH3D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceDm2Sn2GMu3D (TH3D &surface, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceDm2Sn2GMu3DBoth (TH3D &surface, NuMMParameters mmPars)
virtual Double_t TransGMuFixDm2Sn2 (TH1D &surface, TH1D &dm2vals, TH1D &sn2vals, NuMMParameters mmPars)
virtual Double_t FillLikelihoodSurfaceLED (TH2D &surface, NuMMParameters mmPars)
virtual void BatchModeOn ()
 Supresses much of the usual output.
virtual void QuietModeOn ()
virtual void QuietModeOff ()
virtual void FillCPT (TH2D &surface)
 Fills a 2D histogram that has dimensions Dm2 (Y) and DM2Bar(X). Uses the center of the passed bins.
NuMMParameters FillAlphaBeta (TH2D &J, TH2D &K, const std::string outfilename)

Protected Member Functions

virtual Double_t FillLikelihoodSurfaceGeneric (TH2D &surface, NuMMParameters mmPars, int Ix=1, int Iy=1, int Nx=0, int Ny=0, TH2D *mmPars_surfaces[]=0, TH2D **converge_surface=0, int par_index=-1, int n_iter=1, double *iter=0)
virtual Double_t FillLikelihoodSurfaceGeneric (TH3D &surface, NuMMParameters mmPars)
NuMMParameters CoarseGridSearch (NuMMParameters mmGrid)
NuMMParameters GridSearch (TH2D &likesurf, NuMMParameters mmGrid, NuMMParameters *Ref=0)
NuMMParameters FineGridSearch (NuMMParameters mmStart)
NuMMParameters CoarseGridSearchNu (NuMMParameters mmGrid)
NuMMParameters GridSearchNu (TH2D &likesurf, NuMMParameters mmGrid, NuMMParameters *Ref=0)
NuMMParameters FineGridSearchNu (NuMMParameters mmStart)
NuMMParameters CoarseGridSearchNSI (NuMMParameters mmGrid)
NuMMParameters GridSearchNSI (TH2D &likesurf, NuMMParameters mmGrid, NuMMParameters *Ref=0)
NuMMParameters FineGridSearchNSI (NuMMParameters mmStart)
NuMMParameters Sn2GridSearch (Double_t Dm2, Double_t Dm2Bar, Int_t resolution=0)
 Does a grid search over sn2/sn2bar, with resolution bins per parameter.

Protected Attributes

Int_t fQuietMode
Double_t fUp
Int_t fNumContourPoints
std::vector< Double_t > fsn2PointsForNeutrinoContour
std::vector< Double_t > fsn2PointsForNeutrinoContourPhysical
std::vector< Double_t > fsn2PointsForNeutrinoContourUnphysical
std::vector< NuMMRun * > * fRuns

Detailed Description

Definition at line 19 of file NuSystFitter.h.


Constructor & Destructor Documentation

NuSystFitter::NuSystFitter (  ) 

Definition at line 34 of file NuSystFitter.cxx.

00035   : fQuietMode(0),
00036     fUp(1.0),
00037     fNumContourPoints(80),
00038     fsn2PointsForNeutrinoContour(),
00039     fsn2PointsForNeutrinoContourPhysical(),
00040     fsn2PointsForNeutrinoContourUnphysical()
00041 {
00042   fRuns = new vector<NuMMRun*>();
00043 
00044   for (Double_t sn2 = 1.0; sn2 > 0.871; sn2 -= 0.02){
00045     fsn2PointsForNeutrinoContourPhysical.push_back(sn2);
00046   }
00047   fsn2PointsForNeutrinoContourPhysical.push_back(0.87);
00048   fsn2PointsForNeutrinoContourPhysical.push_back(0.86);
00049   for (Double_t sn2 = 0.8595; sn2 > 0.7999; sn2 -= 0.0005){
00050     fsn2PointsForNeutrinoContourPhysical.push_back(sn2);
00051   }
00052 
00053   for (Double_t sn2 = 1.0; sn2 < 3; sn2 += 0.02){
00054     fsn2PointsForNeutrinoContourUnphysical.push_back(sn2);
00055   }
00056 
00057   fsn2PointsForNeutrinoContour = fsn2PointsForNeutrinoContourPhysical;
00058 }

NuSystFitter::~NuSystFitter (  )  [virtual]

Definition at line 61 of file NuSystFitter.cxx.

00062 {
00063 }


Member Function Documentation

void NuSystFitter::BatchModeOn (  )  [virtual]

Supresses much of the usual output.

Definition at line 677 of file NuSystFitter.cxx.

References fQuietMode, and QuietModeOn().

Referenced by FeldmanSterile::Calculate2DdeltaChi2(), NuFCFitter::NuFCFitter(), NuFCFitterNSI::NuFCFitterNSI(), NuFCFitterNSINu::NuFCFitterNSINu(), and NuFCFitterNSINubar::NuFCFitterNSINubar().

00678 {
00679   QuietModeOn();
00680   fQuietMode = 2;
00681 }

const TGraph * NuSystFitter::Chi2SliceDm2 ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minval = 1e-3,
const Double_t  maxval = 3e-3 
) [virtual]

Definition at line 3192 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), fQuietMode, minimum, and NuMMParameters::VectorParameters().

03193 {
03194         Double_t increment = (maxval - minval) / (Double_t)numpoints;
03195 
03196 
03197         if (!fQuietMode) cout << "Generating Likelyhood surface for dm2" << endl;
03198 
03199         // Set the parameters object
03200         NuMMParameters mmPars(mmParsC);
03201 
03202         // Make the TGraph to fill
03203         TGraph *chisurf = new TGraph(numpoints);
03204         chisurf->SetNameTitle("Chi2Dm2", "#Delta#chi^{2} Surface for #Deltam^{2}");
03205         chisurf->GetXaxis()->SetTitle("#Deltam^{2}");
03206 
03207         Double_t dm2value = 0, like = 0, minimum = 0;
03208 
03209         // We want to scan over the sin2bar surface and generate a likelyhood plot
03210         for (Int_t i = 0; i < numpoints; i++)
03211         {
03212                 dm2value = minval + increment*i;
03213                 mmPars.Dm2(dm2value);
03214 
03215                 // Generate the likelyhood....
03216                 like = (*this)(mmPars.VectorParameters());
03217 
03218                 // Fill the array
03219 //              dm2s[i] = dm2value;
03220 //              likelyhoods[i] = like;
03221 
03222                 chisurf->SetPoint(i, dm2value, like);
03223 
03224                 // Is this the lowest value?
03225                 if (i == 0)
03226                         minimum = like;
03227                 else
03228                         minimum = like < minimum ? like : minimum;
03229         }
03230 
03231         // Offset the data points
03232         for (Int_t i = 0; i < numpoints; i++)
03233         {
03234                 Double_t x, y;
03235                 chisurf->GetPoint(i, x, y);
03236                 chisurf->SetPoint(i, x, y-minimum);
03237         }
03238 
03239         // Now
03240         if (!fQuietMode) cout << "Surface Generation Done." << endl;
03241 
03242         // Return the completed TGraph
03243         return chisurf;
03244 }

const TGraph * NuSystFitter::Chi2SliceDm2bar ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minval = 1e-3,
const Double_t  maxval = 3e-3 
) [virtual]

Definition at line 760 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), fQuietMode, minimum, and NuMMParameters::VectorParameters().

00761 {
00762         Double_t increment = (maxval - minval) / (Double_t)numpoints;
00763 
00764 
00765         if (!fQuietMode) cout << "Generating Likelyhood surface for dm2bar" << endl;
00766 
00767         // Set the parameters object
00768         NuMMParameters mmPars(mmParsC);
00769 
00770         // Make the TGraph to fill
00771         TGraph *chisurf = new TGraph(numpoints);
00772         chisurf->SetNameTitle("Chi2Dm2Bar", "#Delta#chi^{2} Surface for #Delta#barm^{2}");
00773         chisurf->GetXaxis()->SetTitle("#Delta#barm^{2}");
00774 
00775         Double_t dm2value = 0, like = 0, minimum = 0;
00776 
00777         // We want to scan over the sin2bar surface and generate a likelyhood plot
00778         for (Int_t i = 0; i < numpoints; i++)
00779         {
00780                 dm2value = minval + increment*i;
00781                 mmPars.Dm2Bar(dm2value);
00782 
00783                 // Generate the likelyhood....
00784                 like = (*this)(mmPars.VectorParameters());
00785 
00786                 // Fill the array
00787 //              dm2s[i] = dm2value;
00788 //              likelyhoods[i] = like;
00789 
00790                 chisurf->SetPoint(i, dm2value, like);
00791 
00792                 // Is this the lowest value?
00793                 if (i == 0)
00794                         minimum = like;
00795                 else
00796                         minimum = like < minimum ? like : minimum;
00797         }
00798 
00799         // Offset the data points
00800         for (Int_t i = 0; i < numpoints; i++)
00801         {
00802                 Double_t x, y;
00803                 chisurf->GetPoint(i, x, y);
00804                 chisurf->SetPoint(i, x, y-minimum);
00805         }
00806 
00807         // Now
00808         if (!fQuietMode) cout << "Surface Generation Done." << endl;
00809 
00810         // Return the completed TGraph
00811         return chisurf;
00812 }

const TGraph * NuSystFitter::Chi2SliceSin2bar ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minval = 0.8,
const Double_t  maxval = 1.2 
) [virtual]

Definition at line 706 of file NuSystFitter.cxx.

References fQuietMode, minimum, NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

00707 {
00708         Double_t increment = (maxval - minval) / (Double_t)numpoints;
00709 
00710         if (!fQuietMode) cout << "Generating Likelyhood surface for sin2bar" << endl;
00711 
00712         // Set the parameters object
00713         NuMMParameters mmPars(mmParsC);
00714 
00715         Double_t sn2value = 0, like = 0, minimum = 0;
00716 
00717         // Make the TGraph to fill
00718         TGraph *chisurf = new TGraph(numpoints);
00719         chisurf->SetNameTitle("Chi2Sin2Bar", "#Delta#chi^{2} Surface for Sin^{2}2#bar#theta");
00720         chisurf->GetXaxis()->SetTitle("Sin^{2}2#bar#theta");
00721 
00722         // We want to scan over the sin2bar surface and generate a likelyhood plot
00723         for (Int_t i = 0; i < numpoints; i++)
00724         {
00725                 sn2value = minval + increment*i;
00726                 mmPars.Sn2Bar(sn2value);
00727 
00728                 // Generate the likelyhood....
00729                 like = (*this)(mmPars.VectorParameters());
00730 
00731                 // Fill the array
00732                 //sin2s[i] = sn2value;
00733                 //likelyhoods[i] = like;
00734                 chisurf->SetPoint(i, sn2value, like);
00735 
00736                 // Is this the lowest value?
00737                 if (i == 0)
00738                         minimum = like;
00739                 else
00740                         minimum = like < minimum ? like : minimum;
00741         }
00742 
00743         // Offset the data points
00744         for (Int_t i = 0; i < numpoints; i++)
00745         {
00746                 Double_t x, y;
00747                 chisurf->GetPoint(i, x, y);
00748                 chisurf->SetPoint(i, x, y-minimum);
00749         }
00750 
00751         // Now
00752         if (!fQuietMode) cout << "Surface Generation Done." << endl;
00753 
00754         // Return the completed TGraph
00755         return chisurf;
00756 }

const TGraph * NuSystFitter::Chi2ValleyBar ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minSin2bar = 0.,
const Double_t  maxSin2bar = 2. 
) [virtual]

Definition at line 815 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixNCBackground(), NuMMParameters::FixNorm(), NuMMParameters::FixShwScale(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), fQuietMode, Minimise(), minimum, NuMMParameters::MinuitPass(), NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

00816 {
00817     NuMMParameters mmPars(mmParsC);
00818 
00819     //  mmPars.Dm2Bar(2.5e-3);
00820     // mmPars.Sn2Bar(1.0);
00821     mmPars.FixNorm();
00822     mmPars.FixShwScale();
00823     mmPars.FixNCBackground();
00824     mmPars.FixDm2();
00825     mmPars.FixSn2();
00826     mmPars.FixDm2Bar(false);
00827     mmPars.FixSn2Bar();
00828 
00829     //  Double_t sin2[points], likelyhood[points], mass[points];
00830     Double_t sinval=0, likeval=0, minimum=0;
00831     Double_t interval =  (max - min) / (points -1);;
00832 
00833     if(!fQuietMode) cout << "Generating Valley trace plot..." << endl;
00834 
00835     TGraph *sin2like = new TGraph(points);
00836 
00837     for (Int_t i = 0; i < points; i++) {
00838 
00839         //        cout << endl << "Sin2 " << i+1 << " / " << points << endl;
00840         // Set up the parameter structure
00841         sinval = min + interval*i;
00842         mmPars.Sn2Bar(sinval);
00843 
00844         NuMMParameters mmFit;
00845         if (sinval == 0) {
00846             mmFit = mmPars;
00847         }
00848         else {
00849             // Now, minimize the likelyhood for this value
00850             mmFit = this->Minimise(mmPars);
00851             int failCount = 0;
00852             while (!mmFit.MinuitPass()) {
00853                 if (failCount > 10) break;
00854                 failCount++;
00855                 mmPars.Dm2Bar(mmPars.Dm2Bar()*1.1);
00856                 cout << "Starting at dm2= " << mmPars.Dm2Bar() << endl;
00857                 mmFit = this->Minimise(mmPars);
00858             }
00859             mmPars.Dm2Bar(mmParsC.Dm2Bar());
00860             if (failCount > 10) continue;
00861         }
00862 
00863         likeval = (*this)(mmFit.VectorParameters());
00864         sin2like->SetPoint(i, sinval, likeval);
00865 
00866         // Calculate Minimum
00867         if (i == 0)
00868             minimum = likeval;
00869         else
00870             minimum = likeval < minimum ? likeval : minimum;
00871         //              cout << "Minimum is now " << minimum << endl;
00872     }
00873 
00874     if (!fQuietMode) cout << endl << "Adjusting Graph Minimum" << endl;
00875     // Offset the data points
00876     for (Int_t i = 0; i < points; i++)
00877         {
00878         Double_t x, y;
00879         sin2like->GetPoint(i, x, y);
00880         //              cout << "From: " << x << ", " << y << " to ";
00881         sin2like->SetPoint(i, x, y-minimum);
00882         //              cout << x << ", " << y-minimum << endl;
00883         }
00884 
00885     if (!fQuietMode) {
00886         cout << endl << "Valley trace plot generation completed" << endl;
00887         cout << "   Minimum found was " << minimum << endl;
00888     }
00889 
00890     // Make the graphs
00891     return sin2like;
00892 }

const TGraph * NuSystFitter::Chi2ValleyBarDm2Bar ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minDm2bar = 0.,
const Double_t  maxDm2bar = 10e-3 
) [virtual]

Definition at line 895 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainPhysical(), NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixNCBackground(), NuMMParameters::FixNorm(), NuMMParameters::FixShwScale(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), fQuietMode, Minimise(), minimum, NuMMParameters::MinuitPass(), NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

00899 {
00900     NuMMParameters mmPars(mmParsC);
00901 
00902     mmPars.FixNorm();
00903     mmPars.FixShwScale();
00904     mmPars.FixNCBackground();
00905     mmPars.FixDm2();
00906     mmPars.FixSn2();
00907     mmPars.FixDm2Bar();
00908     mmPars.FixSn2Bar(false);
00909     mmPars.ConstrainPhysical();
00910 
00911 
00912     //  Double_t sin2[points], likelyhood[points], mass[points];
00913     Double_t dmval=0, likeval=0, minimum=0;
00914     Double_t interval =  (max - min) / (points -1);;
00915 
00916     if(!fQuietMode) cout << "Generating Valley trace plot for dm2..." << endl;
00917 
00918     TGraph *dm2like = new TGraph(points);
00919     dm2like->SetName("valleyplot_mass");
00920 
00921     for (Int_t i = 0; i < points; i++) {
00922         if (i % 10 == 0) {
00923             cout << i << "/" << points << endl;
00924         }
00925         // Set up the parameter structure
00926         dmval = min + interval*i;
00927         mmPars.Dm2Bar(dmval);
00928 
00929         NuMMParameters mmFit;
00930         if (dmval == 0) {
00931             mmFit = mmPars;
00932         }
00933         else {
00934             // Now, minimize the likelyhood for this value
00935             mmFit = this->Minimise(mmPars);
00936             int failCount = 0;
00937             while (!mmFit.MinuitPass()) {
00938                 if (failCount > 10) break;
00939                 failCount++;
00940                 mmPars.Sn2Bar(mmPars.Sn2Bar()*0.9);
00941                 cout << "Starting at sn2= " << mmPars.Sn2Bar() << endl;
00942                 mmFit = this->Minimise(mmPars);
00943             }
00944             mmPars.Sn2Bar(mmParsC.Sn2Bar());
00945             if (failCount > 10) continue;
00946         }
00947         //NuMMParameters mmFit(mmPars);
00948 
00949         likeval = (*this)(mmFit.VectorParameters());
00950         dm2like->SetPoint(i, dmval, likeval);
00951 
00952         // Calculate Minimum
00953         if (i == 0)
00954             minimum = likeval;
00955         else
00956             minimum = likeval < minimum ? likeval : minimum;
00957     }
00958 
00959     if (!fQuietMode) cout << endl << "Adjusting Graph Minimum" << endl;
00960 
00961     // Offset the data points
00962     for (Int_t i = 0; i < points; i++) {
00963         Double_t x, y;
00964         dm2like->GetPoint(i, x, y);
00965         dm2like->SetPoint(i, x, y-minimum);
00966     }
00967 
00968     if (!fQuietMode) {
00969         cout << endl << "Valley trace plot generation completed" << endl;
00970         cout << "   Minimum found was " << minimum << endl;
00971     }
00972 
00973     return dm2like;
00974 }

const TGraph * NuSystFitter::Chi2ValleyDm2 ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minDm2 = 0.,
const Double_t  maxDm2 = 10e-3 
) [virtual]

Definition at line 978 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainPhysical(), NuMMParameters::Dm2(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixNCBackground(), NuMMParameters::FixNorm(), NuMMParameters::FixShwScale(), NuMMParameters::FixSn2Bar(), fQuietMode, Minimise(), minimum, NuUtilities::ProgressBar(), NuMMParameters::ReleaseSn2(), and NuMMParameters::VectorParameters().

00982 {
00983     NuMMParameters mmPars(mmParsC);
00984 
00985     mmPars.FixNorm();
00986     mmPars.FixShwScale();
00987     mmPars.FixNCBackground();
00988     mmPars.FixDm2();
00989     mmPars.ReleaseSn2();
00990     mmPars.FixDm2Bar();
00991     mmPars.FixSn2Bar();
00992     mmPars.ConstrainPhysical();
00993 
00994 
00995     //  Double_t sin2[points], likelyhood[points], mass[points];
00996     Double_t dmval=0, likeval=0, minimum=0;
00997     Double_t interval =  (max - min) / (points -1);;
00998 
00999     if(!fQuietMode) cout << "Generating Valley trace plot for dm2..." << endl;
01000 
01001     TGraph *dm2like = new TGraph(points);
01002     dm2like->SetName("valleyplot_mass");
01003 
01004     for (Int_t i = 0; i < points; i++)
01005     {
01006         // Print a progress bar
01007         NuUtilities::ProgressBar(i, points-1,5);
01008 
01009         // Set up the parameter structure
01010         dmval = min + interval*i;
01011         mmPars.Dm2(dmval);
01012 
01013         // Now, minimize the likelyhood for this value
01014         NuMMParameters mmFit = this->Minimise(mmPars);
01015         //NuMMParameters mmFit(mmPars);
01016 
01017         likeval = (*this)(mmFit.VectorParameters());
01018         dm2like->SetPoint(i, dmval, likeval);
01019 
01020         // Calculate Minimum
01021         if (i == 0)
01022             minimum = likeval;
01023         else
01024             minimum = likeval < minimum ? likeval : minimum;
01025     }
01026 
01027     if (!fQuietMode) cout << endl << "Adjusting Graph Minimum" << endl;
01028 
01029     // Offset the data points
01030     for (Int_t i = 0; i < points; i++)
01031     {
01032         Double_t x, y;
01033         dm2like->GetPoint(i, x, y);
01034         dm2like->SetPoint(i, x, y-minimum);
01035     }
01036 
01037     if (!fQuietMode) {
01038         cout << endl << "Valley trace plot generation completed" << endl;
01039         cout << "   Minimum found was " << minimum << endl;
01040     }
01041 
01042     return dm2like;
01043 }

const TGraph * NuSystFitter::Chi2ValleySn2 ( const NuMMParameters mmPars,
const Int_t  numpoints = 100,
const Double_t  minSn2 = 0.,
const Double_t  maxSn2 = 1. 
) [virtual]

Definition at line 1047 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainPhysical(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixNCBackground(), NuMMParameters::FixNorm(), NuMMParameters::FixShwScale(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), fQuietMode, Minimise(), minimum, NuUtilities::ProgressBar(), NuMMParameters::ReleaseDm2(), NuMMParameters::Sn2(), and NuMMParameters::VectorParameters().

01051 {
01052     NuMMParameters mmPars(mmParsC);
01053 
01054     mmPars.FixNorm();
01055     mmPars.FixShwScale();
01056     mmPars.FixNCBackground();
01057     mmPars.ReleaseDm2();
01058     mmPars.FixSn2();
01059     mmPars.FixDm2Bar();
01060     mmPars.FixSn2Bar();
01061     mmPars.ConstrainPhysical();
01062 
01063 
01064     //  Double_t sin2[points], likelyhood[points], mass[points];
01065     Double_t snval=0, likeval=0, minimum=0;
01066     Double_t interval =  (max - min) / (points -1);;
01067 
01068     if(!fQuietMode) cout << "Generating Valley trace plot for sn2..." << endl;
01069 
01070     TGraph *sn2like = new TGraph(points);
01071     sn2like->SetName("valleyplot_sin");
01072 
01073     for (Int_t i = 0; i < points; i++)
01074     {
01075         // Print a progress bar
01076         NuUtilities::ProgressBar(i, points-1,5);
01077 
01078         // Set up the parameter structure
01079         snval = min + interval*i;
01080         mmPars.Sn2(snval);
01081 
01082         // Now, minimize the likelyhood for this value
01083         NuMMParameters mmFit = this->Minimise(mmPars);
01084         //NuMMParameters mmFit(mmPars);
01085 
01086         likeval = (*this)(mmFit.VectorParameters());
01087         sn2like->SetPoint(i, snval, likeval);
01088 
01089         // Calculate Minimum
01090         if (i == 0)
01091             minimum = likeval;
01092         else
01093             minimum = likeval < minimum ? likeval : minimum;
01094     }
01095 
01096     if (!fQuietMode) cout << endl << "Adjusting Graph Minimum" << endl;
01097 
01098     // Offset the data points
01099     for (Int_t i = 0; i < points; i++)
01100     {
01101         Double_t x, y;
01102         sn2like->GetPoint(i, x, y);
01103         sn2like->SetPoint(i, x, y-minimum);
01104     }
01105 
01106     if (!fQuietMode) {
01107         cout << endl << "Valley trace plot generation completed" << endl;
01108         cout << "   Minimum found was " << minimum << endl;
01109     }
01110 
01111     return sn2like;
01112 }

NuMMParameters NuSystFitter::CoarseGridSearch ( NuMMParameters  mmGrid  )  [protected]

Definition at line 2606 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), MuELoss::e, GridSearch(), Msg::kInfo, MSG, and NuMMParameters::Sn2Bar().

Referenced by FCMinimise().

02607 {
02608   // The width and height of the grid
02609   const Int_t gridWidth = 5;
02610   const Int_t gridHeight = 15;
02611 
02612   // Generate the surface
02613   TH2D likesurf("likesurf", "likesurf", gridWidth, 0.2, 1, gridHeight, 2e-3, 30e-3);
02614   // TH2D likesurf("likesurf", "likesurf", gridWidth, 0.2, 1, 998*10, 2e-3, 1000e-3);
02615   mmGrid = GridSearch(likesurf, mmGrid);
02616   // cout << "Found dm2: " << mmGrid.Dm2Bar() << ", sn2: " << mmGrid.Sn2Bar() << endl;
02617 
02618   // Now do a log search
02619   vector<Double_t> bins;
02620   for (Double_t i = -1.5; i <= -0.; i += 0.1) {
02621     bins.push_back(pow(10, i));
02622   }
02623   Int_t bincount = bins.size()-1;
02624   TH2D logsurf("logsurf", "logsurf", gridWidth, 0.2, 1, bincount, &(bins.front()));
02625   mmGrid = GridSearch(logsurf, mmGrid, &mmGrid);
02626 
02627   // This was code to do a much sparser grid search instead of the log search
02628   // TH2D sparser("sparser", "sparser", gridWidth, 0.2, 1, 400, 30e-3, 1000e-3);
02629   //  mmGrid = GridSearch(sparser, mmGrid, &mmGrid);
02630 
02631   // Output the minimum
02632   MSG("NuSystFitter", Msg::kInfo) << "Coarse grid search found sn2: "
02633     << mmGrid.Sn2Bar() << ", dm: " << mmGrid.Dm2Bar() << '\n';
02634 
02635   return mmGrid;
02636 }

NuMMParameters NuSystFitter::CoarseGridSearchNSI ( NuMMParameters  mmGrid  )  [protected]

Definition at line 2931 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), MuELoss::e, NuMMParameters::Epsilon(), GridSearchNSI(), Msg::kInfo, and MSG.

Referenced by FCMinimiseNSI().

02932 {
02933   // The width and height of the grid
02934   const Int_t gridWidth = 5;
02935   const Int_t gridHeight = 15;
02936 
02937   // Generate the surface
02938   TH2D likesurf("likesurf", "likesurf", gridWidth, 0.01, 0.5, gridHeight, 0e-3, 10e-3);
02939   mmGrid = GridSearchNSI(likesurf, mmGrid);
02940   cout << "Found dm2: " << mmGrid.Dm2() << ", epsilon: " << mmGrid.Epsilon() << endl;
02941 
02942   // Now do a log search
02943   vector<Double_t> bins;
02944   for (Double_t i = -1.5; i <= -0.; i += 0.1) {
02945     bins.push_back(pow(10, i));
02946   }
02947   Int_t bincount = bins.size()-1;
02948   TH2D logsurf("logsurf", "logsurf", gridWidth, 0.01, 0.5, bincount, &(bins.front()));
02949   mmGrid = GridSearchNSI(logsurf, mmGrid, &mmGrid);
02950 
02951   // Output the minimum
02952   MSG("NuSystFitter", Msg::kInfo) << "Coarse grid search found epsilon: "
02953     << mmGrid.Epsilon() << ", dm: " << mmGrid.Dm2() << '\n';
02954 
02955   return mmGrid;
02956 }

NuMMParameters NuSystFitter::CoarseGridSearchNu ( NuMMParameters  mmGrid  )  [protected]

Definition at line 2766 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), MuELoss::e, GridSearchNu(), Msg::kInfo, MSG, and NuMMParameters::Sn2().

Referenced by FCMinimiseNu().

02767 {
02768   // The width and height of the grid
02769   const Int_t gridWidth = 5;
02770   const Int_t gridHeight = 15;
02771   
02772   // Generate the surface
02773   TH2D likesurf("likesurf", "likesurf", gridWidth, 0.2, 1, gridHeight, 2e-3, 30e-3);
02774   // TH2D likesurf("likesurf", "likesurf", gridWidth, 0.2, 1, 998*10, 2e-3, 1000e-3);
02775   mmGrid = GridSearchNu(likesurf, mmGrid);
02776   // cout << "Found dm2: " << mmGrid.Dm2Bar() << ", sn2: " << mmGrid.Sn2Bar() << endl;
02777   
02778   // Now do a log search
02779   vector<Double_t> bins;
02780   for (Double_t i = -1.5; i <= -0.; i += 0.1) {
02781     bins.push_back(pow(10, i));
02782   }
02783   Int_t bincount = bins.size()-1;
02784   TH2D logsurf("logsurf", "logsurf", gridWidth, 0.2, 1, bincount, &(bins.front()));
02785   mmGrid = GridSearchNu(logsurf, mmGrid, &mmGrid);
02786   
02787   // This was code to do a much sparser grid search instead of the log search
02788   // TH2D sparser("sparser", "sparser", gridWidth, 0.2, 1, 400, 30e-3, 1000e-3);
02789   //  mmGrid = GridSearch(sparser, mmGrid, &mmGrid);
02790   
02791   // Output the minimum
02792   MSG("NuSystFitter", Msg::kInfo) << "Coarse grid search found sn2: "
02793   << mmGrid.Sn2() << ", dm: " << mmGrid.Dm2() << endl;
02794   
02795   return mmGrid;
02796 }

const std::pair< Double_t, Double_t > NuSystFitter::DmScanForContour ( const NuMMParameters mmPars,
Double_t  targetChi2 
) [virtual]

Definition at line 384 of file NuSystFitter.cxx.

References NuMMParameters::AreAllParametersFixed(), NuMMParameters::Dm2(), MuELoss::e, NuMMParameters::FixDm2(), NuMMParameters::NCBackgroundScale(), NuMMParameters::Normalisation(), NuMMParameters::ShwEnScale(), NuMMParameters::Sn2(), and NuMMParameters::VectorParameters().

00386 {
00387   std::pair<Double_t,Double_t> errors;
00388 
00389   NuMMParameters bestFitPoint = this->Minimise(mmPars);
00390   Double_t bestFitChi2 = (*this)(bestFitPoint.VectorParameters());
00391   Double_t bestFitDm2 = bestFitPoint.Dm2();
00392 
00393   if (bestFitChi2 > targetChi2){
00394     errors.first = -1.0;
00395     errors.second = -1.0;
00396     return errors;
00397   }
00398 
00399   cout << "Best chi2: " << bestFitChi2
00400        << " at dm2 = " << bestFitDm2
00401        << ", sn2 = " << bestFitPoint.Sn2()
00402        << ", norm = " << bestFitPoint.Normalisation()
00403        << ", NC = " << bestFitPoint.NCBackgroundScale()
00404        << ", shwEn = " << bestFitPoint.ShwEnScale()
00405        << endl
00406        << ". Target chi2: " << targetChi2
00407        << endl;
00408 
00409   for (Int_t direction = -1; direction < 2; direction += 2){
00410 
00411     cout << "Direction = " << direction << endl;
00412 
00413     NuMMParameters scanMMPars = mmPars;
00414     scanMMPars.Dm2(bestFitDm2);
00415     scanMMPars.FixDm2();
00416 
00417     Double_t incrementArray[9] = {0.1e-3,
00418                                   0.03e-3,
00419                                   0.01e-3,
00420                                   0.003e-3,
00421                                   0.001e-3,
00422                                   0.0003e-3,
00423                                   0.0001e-3,
00424                                   0.00003e-3,
00425                                   0.00001e-3};
00426     Double_t firstChi2 = bestFitChi2;
00427     Double_t firstDm2 = bestFitDm2;
00428     Double_t currentChi2 = -1.0;
00429     Double_t lastChi2 = -1.0;
00430     Double_t lastDm2 = -1.0;
00431     for (Int_t incCount = 0; incCount<9; ++incCount){
00432       Double_t increment = incrementArray[incCount];
00433 
00434       currentChi2 = firstChi2;
00435       scanMMPars.Dm2(firstDm2);
00436 
00437       cout << "Starting while with firstChi2 = " << firstChi2
00438            << ", firstDm2 = " << firstDm2
00439            << ", increment = " << increment
00440            << endl;
00441 
00442       while (currentChi2 < targetChi2){
00443         lastChi2 = currentChi2;
00444         lastDm2 = scanMMPars.Dm2();
00445 
00446         scanMMPars.Dm2(scanMMPars.Dm2() + direction*increment);
00447         cout << "Trying dm2 = " << scanMMPars.Dm2() << "... ";
00448         flush(cout);
00449         if (scanMMPars.AreAllParametersFixed()){
00450           currentChi2 = (*this)(scanMMPars.VectorParameters());
00451         }
00452         else{
00453           currentChi2 = (*this)(this->Minimise(scanMMPars).VectorParameters());
00454         }
00455         cout << "chi2 = " << currentChi2 << endl;
00456       }
00457 
00458       firstDm2 = lastDm2;
00459       firstChi2 = lastChi2;
00460 
00461     }
00462 
00463     Double_t chi2Difference = lastChi2 - currentChi2;
00464     Double_t dm2Difference = lastDm2 - scanMMPars.Dm2();
00465     Double_t value = -1.0;
00466     cout << "chi2Difference = " << chi2Difference << endl;
00467     cout << "dm2Difference = " << dm2Difference << endl;
00468     if ((fabs(chi2Difference)<1e-12) || (fabs(dm2Difference)<1e-12)){
00469       cout << "Chi2 or Dm2 difference is zero. Supplying currentChi2."
00470            << endl;
00471       value = scanMMPars.Dm2();
00472     }
00473     else{
00474       Double_t fractionBetween = (targetChi2 - currentChi2)/chi2Difference;
00475       value = fractionBetween*dm2Difference + scanMMPars.Dm2();
00476     }
00477     cout << "Error is between " << lastDm2
00478          << " and " << scanMMPars.Dm2()
00479          << ", with chi2 between " << lastChi2
00480          << " and " << currentChi2
00481          << endl;
00482     cout << "Value is " << value << endl;
00483     if (-1 == direction){errors.first = value;}
00484     if (1 == direction){errors.second = value;}
00485   }
00486 
00487   cout << "Final answer is dm2 = " << bestFitDm2
00488        << " + " << errors.second - bestFitDm2
00489        << " - " << errors.first - bestFitDm2
00490        << endl;
00491   return errors;
00492 }

const NuMMParameters NuSystFitter::FCMinimise ( const NuMMParameters mmPars  )  [virtual]

More robustly minimizes by doing a wide grid search first - same as the FC fitting.

Definition at line 2572 of file NuSystFitter.cxx.

References CoarseGridSearch(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), FineGridSearch(), Msg::kInfo, Minimise(), NuMMParameters::MinuitPass(), MSG, NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

Referenced by NuFCFitterNuMuBar::Fit().

02573 {
02574   // Start with a coarse grid search, to seed minuit
02575   const NuMMParameters mmStart = CoarseGridSearch(mmPars);
02576 
02577   // Now run the normal minuit procedure with this as the starting point
02578   NuMMParameters mmFit = Minimise(mmStart);
02579 
02580   // If minuit failed, do a finer grid search around the coarse starting point
02581   // Did minuit succeed?
02582   // if (mmFit.MinuitPass())
02583   // {
02584 
02585   // } else {
02586   if (!mmFit.MinuitPass())
02587   {
02588       // Minuit failed to fit. Use a more intensive grid search,
02589       // starting on the minuit starting point (that we get from
02590       // a grid search)
02591       MSG("NuSystFitter", Msg::kInfo) << "Minuit fit failed. Doing Finer Grid search" << endl;
02592       mmFit = FineGridSearch(mmStart);
02593   } else {
02594       // Tell them what minuit found
02595       MSG("NuSystFitter", Msg::kInfo)
02596         << "Minuit found: dm2: " << mmFit.Dm2()
02597         << " sn2: " << mmFit.Sn2()
02598         << " dm2bar: " << mmFit.Dm2Bar()
02599         << " sn2bar: " << mmFit.Sn2Bar()
02600         << '\n';
02601   }
02602   return mmFit;
02603 }

const NuMMParameters NuSystFitter::FCMinimiseNSI ( const NuMMParameters mmPars  )  [virtual]

Definition at line 2891 of file NuSystFitter.cxx.

References CoarseGridSearchNSI(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Epsilon(), FineGridSearchNSI(), Msg::kInfo, Minimise(), NuMMParameters::MinuitPass(), MSG, NuMMParameters::PrintStatus(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

02892 {
02893   // Start with a coarse grid search, to seed minuit
02894   mmPars.PrintStatus();
02895   const NuMMParameters mmStart = CoarseGridSearchNSI(mmPars);
02896 
02897   // Now run the normal minuit procedure with this as the starting point
02898   mmStart.PrintStatus();
02899   NuMMParameters mmFit = Minimise(mmStart);
02900 
02901   // Test seeding directly with no coarse search
02902   //NuMMParameters mmFit = Minimise(mmPars);
02903 
02904   // If minuit failed, do a finer grid search around the coarse starting point
02905   // Did minuit succeed?
02906   // if (mmFit.MinuitPass())
02907   // {
02908 
02909   // } else {
02910   if (!mmFit.MinuitPass())
02911   {
02912       // Minuit failed to fit. Use a more intensive grid search,
02913       // starting on the minuit starting point (that we get from
02914       // a grid search)
02915       MSG("NuSystFitter", Msg::kInfo) << "Minuit fit failed. Doing Finer Grid search" << endl;
02916       mmFit = FineGridSearchNSI(mmStart);
02917   } else {
02918       // Tell them what minuit found
02919       MSG("NuSystFitter", Msg::kInfo)
02920         << "Minuit found: dm2: " << mmFit.Dm2()
02921         << " sn2: " << mmFit.Sn2()
02922         << " dm2bar: " << mmFit.Dm2Bar()
02923         << " sn2bar: " << mmFit.Sn2Bar()
02924         << "epsilon: " << mmFit.Epsilon()
02925         << '\n';
02926   }
02927   return mmFit;
02928 }

const NuMMParameters NuSystFitter::FCMinimiseNu ( const NuMMParameters mmPars  )  [virtual]

Definition at line 2732 of file NuSystFitter.cxx.

References CoarseGridSearchNu(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), FineGridSearchNu(), Msg::kInfo, Minimise(), NuMMParameters::MinuitPass(), MSG, NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

02733 {
02734   // Start with a coarse grid search, to seed minuit
02735   const NuMMParameters mmStart = CoarseGridSearchNu(mmPars);    
02736   
02737   // Now run the normal minuit procedure with this as the starting point
02738   NuMMParameters mmFit = Minimise(mmStart);
02739   
02740   // If minuit failed, do a finer grid search around the coarse starting point
02741   // Did minuit succeed?
02742   // if (mmFit.MinuitPass())
02743   // {
02744   
02745   // } else {
02746   if (!mmFit.MinuitPass())
02747   {
02748     // Minuit failed to fit. Use a more intensive grid search,
02749     // starting on the minuit starting point (that we get from
02750     // a grid search)
02751     MSG("NuSystFitter", Msg::kInfo) << "Minuit fit failed. Doing Finer Grid search" << endl;
02752     mmFit = FineGridSearchNu(mmStart);
02753   } else {
02754     // Tell them what minuit found
02755     MSG("NuSystFitter", Msg::kInfo)
02756     << "Minuit found: dm2: " << mmFit.Dm2()
02757     << " sn2: " << mmFit.Sn2()
02758     << " dm2bar: " << mmFit.Dm2Bar()
02759     << " sn2bar: " << mmFit.Sn2Bar()
02760     << '\n';
02761   }
02762   return mmFit;
02763 }

NuMMParameters NuSystFitter::FillAlphaBeta ( TH2D &  J,
TH2D &  K,
const std::string  outfilename 
)

Definition at line 3386 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), minimum, NuUtilities::ProgressBar(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

03387 {
03388   // characteristics of each surface
03389   Int_t NpointsJX = J.GetNbinsX();
03390   Int_t NpointsJY = J.GetNbinsY();
03391   Int_t NpointsKX = K.GetNbinsX();
03392   Int_t NpointsKY = K.GetNbinsY();
03393 
03394   Double_t widthX = J.GetXaxis()->GetBinWidth(1);
03395   Double_t widthY = J.GetYaxis()->GetBinWidth(1);
03396 
03397   Double_t minKX = K.GetXaxis()->GetBinCenter(1);
03398   Double_t maxKX = K.GetXaxis()->GetBinCenter(NpointsKX);
03399   Double_t minKY = K.GetYaxis()->GetBinCenter(1);
03400   Double_t maxKY = K.GetYaxis()->GetBinCenter(NpointsKY);
03401 
03402   Double_t minJX = J.GetXaxis()->GetBinCenter(1);
03403   Double_t maxJX = J.GetXaxis()->GetBinCenter(NpointsJX);
03404   Double_t minJY = J.GetYaxis()->GetBinCenter(1);
03405   Double_t maxJY = J.GetYaxis()->GetBinCenter(NpointsJY);
03406 
03407   Double_t minABX = minKX - maxJX + widthX/2;
03408   Double_t maxABX = maxKX - minJX - widthX/2.;
03409   Double_t minABY = minKY - maxJY + widthY/2.;
03410   Double_t maxABY = maxKY - minJY - widthY/2.;
03411 
03412   Int_t NpointsAlpha = (int)((maxABY-minABY)/widthY) + 1;
03413   Int_t NpointsBeta  = (int)((maxABX-minABX)/widthX) + 1;
03414 
03415   // resulting (alpha, beta) surfaces
03416   TH2D *absurf = new TH2D("ab", "(#alpha,#beta) likelihood surface",
03417                           NpointsBeta, minABX, maxABX,
03418                           NpointsAlpha, minABY, maxABY );
03419   absurf->GetXaxis()->SetTitle("#beta");
03420   absurf->GetYaxis()->SetTitle("#alpha");
03421 
03422   TH2D *hDmMin = new TH2D("hDmMin", "Dm2 for Minimum Likehood",
03423                           NpointsBeta, minABX, maxABX,
03424                           NpointsAlpha, minABY, maxABY );
03425   hDmMin->GetXaxis()->SetTitle("#beta");
03426   hDmMin->GetYaxis()->SetTitle("#alpha");
03427 
03428   TH2D *hSnMin = new TH2D("hSnMin", "Sn2 for Minimum Likehood",
03429                           NpointsBeta, minABX, maxABX,
03430                           NpointsAlpha, minABY, maxABY );
03431   hSnMin->GetXaxis()->SetTitle("#beta");
03432   hSnMin->GetYaxis()->SetTitle("#alpha");
03433 
03434   TH2D *hDmMinBar = new TH2D("hDmMinBar", "Dm2Bar for Minimum Likehood",
03435                           NpointsBeta, minABX, maxABX,
03436                           NpointsAlpha, minABY, maxABY );
03437   hDmMinBar->GetXaxis()->SetTitle("#beta");
03438   hDmMinBar->GetYaxis()->SetTitle("#alpha");
03439 
03440   TH2D *hSnMinBar = new TH2D("hSnMinBar", "Sn2Bar for Minimum Likehood",
03441                           NpointsBeta, minABX, maxABX,
03442                           NpointsAlpha, minABY, maxABY );
03443   hSnMinBar->GetXaxis()->SetTitle("#beta");
03444   hSnMinBar->GetYaxis()->SetTitle("#alpha");
03445 
03446   // building the surfaces
03447   Double_t likeTot, minimum;
03448   Double_t dmMin=0, snMin=0, dmMinBar=0, snMinBar=0;
03449 
03450   for (Double_t beta=minABX; beta<=maxABX; beta+=widthX){
03451     Int_t iB = absurf->GetXaxis()->FindFixBin(beta);
03452 
03453     for(Double_t alpha=minABY; alpha<=maxABY; alpha+=widthY){
03454       Int_t iA = absurf->GetYaxis()->FindFixBin(alpha);
03455       minimum = -1;
03456       NuUtilities::ProgressBar(iA*iB, NpointsAlpha*NpointsBeta, 1);
03457 
03458       for(Int_t indexJx=1; indexJx<=NpointsJX; ++indexJx){
03459         for(Int_t indexJy=1; indexJy<=NpointsJY; ++indexJy){
03460 
03461           Double_t xval = J.GetXaxis()->GetBinCenter(indexJx);
03462           Double_t yval = J.GetYaxis()->GetBinCenter(indexJy);
03463 
03464           xval += beta;
03465           if (xval < minKX) continue;
03466           if (xval > maxKX) continue;
03467 
03468           yval += alpha;
03469           if (yval < minKY) continue;
03470           if (yval > maxKY) continue;
03471 
03472           Int_t indexKx = K.GetXaxis()->FindFixBin(xval);
03473           Int_t indexKy = K.GetYaxis()->FindFixBin(yval);
03474           likeTot = J.GetBinContent(indexJx, indexJy) + K.GetBinContent(indexKx, indexKy);
03475 
03476           if ( (minimum==-1) || (likeTot<minimum) ) {
03477             minimum = likeTot;
03478             dmMin    = J.GetYaxis()->GetBinCenter(indexJy);
03479             snMin    = J.GetXaxis()->GetBinCenter(indexJx);
03480             dmMinBar = K.GetYaxis()->GetBinCenter(indexKy);
03481             snMinBar = K.GetXaxis()->GetBinCenter(indexKx);
03482           }
03483         }
03484       }
03485       absurf->SetBinContent(iB, iA, minimum);
03486       hDmMin->SetBinContent(iB, iA, dmMin);
03487       hSnMin->SetBinContent(iB, iA, snMin);
03488       hDmMinBar->SetBinContent(iB, iA, dmMinBar);
03489       hSnMinBar->SetBinContent(iB, iA, snMinBar);
03490     }
03491   }
03492 
03493   // shift by minimum
03494   Int_t minbinX, minbinY, minbinZ;
03495   absurf->GetMinimumBin(minbinX, minbinY, minbinZ);
03496   Double_t minContent = absurf->GetBinContent(minbinX, minbinY);
03497 
03498   for (Int_t iterB=1; iterB<=NpointsBeta; ++iterB){
03499     for (Int_t iterA=1; iterA<=NpointsAlpha; ++iterA){
03500       Double_t content = absurf->GetBinContent(iterB, iterA);
03501       absurf->SetBinContent(iterB, iterA, content-minContent);
03502     }
03503   }
03504 
03505   // min of Alpha Beta surface
03506   Double_t minalpha = absurf->GetYaxis()->GetBinCenter(minbinY);
03507   Double_t minbeta  = absurf->GetXaxis()->GetBinCenter(minbinX);
03508   TGraph *minG = new TGraph(1, &minbeta, &minalpha);
03509 
03510   // save
03511   cout << "Plots saved in " << outfilename.c_str() << endl;
03512   TFile *file = new TFile(outfilename.c_str(), "RECREATE");
03513 
03514   J.Write("J", TObject::kOverwrite);
03515   K.Write("K", TObject::kOverwrite);
03516   absurf->Write("AlphaBeta", TObject::kOverwrite);
03517   hDmMin->Write("dm2Min", TObject::kOverwrite);
03518   hSnMin->Write("sn2Min", TObject::kOverwrite);
03519   hDmMinBar->Write("dm2MinBar", TObject::kOverwrite);
03520   hSnMinBar->Write("sn2MinBar", TObject::kOverwrite);
03521   minG->Write("AlphaBetaMinimum", TObject::kOverwrite);
03522   file->Close(0);
03523 
03524   // done!
03525   NuMMParameters pars;
03526   pars.Dm2( hDmMin->GetBinContent(minbinX, minbinY) );
03527   pars.Sn2( hSnMin->GetBinContent(minbinX, minbinY) );
03528   pars.Dm2Bar( hDmMinBar->GetBinContent(minbinX, minbinY) );
03529   pars.Sn2Bar( hSnMinBar->GetBinContent(minbinX, minbinY) );
03530   return pars;
03531 }

void NuSystFitter::FillCPT ( TH2D &  surface  )  [virtual]

Fills a 2D histogram that has dimensions Dm2 (Y) and DM2Bar(X). Uses the center of the passed bins.

Definition at line 3315 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainPhysical(), MuELoss::e, NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), Msg::kError, Msg::kWarning, Minimise(), minimum, NuMMParameters::MinuitPass(), MSG, NuUtilities::ProgressBar(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), Sn2GridSearch(), and NuMMParameters::VectorParameters().

03316 {
03317   Double_t minimum = 0.0;
03318 
03319   const Double_t bincount = surf.GetNbinsX() * surf.GetNbinsY();
03320 
03321   // Loop over every bin in the histogram
03322   for(Int_t xBin = 1; xBin <= surf.GetNbinsX(); ++xBin) {
03323     for (Int_t yBin = 1; yBin <= surf.GetNbinsY(); ++yBin) {
03324       // Show a progress bar
03325       NuUtilities::ProgressBar(surf.GetNbinsY()*(xBin-1) + yBin, bincount, 5);
03326 
03327       const Double_t dm2 =    surf.GetYaxis()->GetBinCenter(yBin);
03328       const Double_t dm2bar = surf.GetXaxis()->GetBinCenter(xBin);
03329 
03330       // Find the sn2 minimum at this value. First do a grid search
03331       NuMMParameters pars = Sn2GridSearch(dm2, dm2bar);
03332       // Get the likelihood, so we can be certain MINUIT improved the fit
03333       const Double_t likeSparse = (*this)(pars.VectorParameters());
03334 
03335       // Now, minimise from this found point
03336       pars.FixDm2();
03337       pars.FixDm2Bar();
03338       pars.ConstrainPhysical();
03339       NuMMParameters parsFit = Minimise(pars);
03340 
03341       // Handle MINUIT failure
03342       if (!parsFit.MinuitPass()) {
03343         MSG("NuSystFitter",Msg::kWarning)
03344           << "MINUIT fit did not pass when attempting to make Dm2-Dm2Bar plot. " << '\n'
03345           << "Dm2: " << dm2 << ", Dm2Bar: " << dm2bar << '\n'
03346           << "Grid Search returned Sn2: " << pars.Sn2() << ", Sn2Bar: " << pars.Sn2Bar() << '\n'
03347           << "Doing Higher resolution grid search" << endl;
03348         // Just do a higher resolution grid search
03349         parsFit = Sn2GridSearch(dm2, dm2bar, 31);
03350         MSG("NuSystFitter",Msg::kWarning)
03351           << "  Found Sn2: " << parsFit.Sn2() << ", Sn2Bar: " << parsFit.Sn2Bar() << endl;
03352       }
03353 
03354       // We now have a best fit for this point. Get the likelihood
03355       const Double_t like = (*this)(parsFit.VectorParameters());
03356       // Sort out the minimums
03357       if (xBin == 1 && yBin == 1) minimum = like;
03358       if (like < minimum) minimum = like;
03359 
03360       // Validate this found minimum was actually better than the grid search
03361       if (like - likeSparse > 1e-4) {
03362         MSG("NuSystFitter",Msg::kError)
03363           << "Dm2: " << dm2 << ", Dm2Bar: " << dm2bar << '\n'
03364           << "Grid Search  Sn2: " << pars.Sn2() << ", Sn2Bar: " << pars.Sn2Bar() << " = " << likeSparse << '\n'
03365           << "MINUIT Found Sn2: " << parsFit.Sn2() << ", Sn2Bar: " << parsFit.Sn2Bar() << " = " << like << '\n'
03366           << "MINUIT fit was worse than seeding grid search!\n"  << '\n'
03367           << "Since this error is not understood, Skipping this bin." << endl;
03368         continue;
03369       }
03370       // Finally, fill the histogram with this likelihood
03371       surf.SetBinContent(xBin, yBin, like);
03372     }
03373   }
03374 
03375   // Now, shift the whole histogram based on the minimum
03376   for(Int_t xBin = 1; xBin <= surf.GetNbinsX(); ++xBin) {
03377     for (Int_t yBin = 1; yBin <= surf.GetNbinsY(); ++yBin) {
03378       Double_t val = surf.GetBinContent(xBin, yBin);
03379       surf.SetBinContent(xBin, yBin, val-minimum);
03380     }
03381   }
03382 }

Double_t NuSystFitter::FillLikelihoodSurface ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Fills a histogram with likelihood values. This function accepts an arbitrarily binned histogram, and fills the values according to the likelihood. It returns the value of the minimum likelihood. (This does the CC neutirno spectrum).

Definition at line 1385 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNuMus(), and FillLikelihoodSurfaceGeneric().

Referenced by GridSearchNu().

01387 {
01388   mmPars.ContourForNuMus();
01389   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01390 }

Double_t NuSystFitter::FillLikelihoodSurfaceBar ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Fills a histogram with likelihood values. This function accepts an arbitrarily binned histogram, and fills the values according to the likelihood. It returns the value of the minimum likelihood.

Definition at line 1394 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNuBars(), and FillLikelihoodSurfaceGeneric().

Referenced by NuFCFitter::Archive(), GridSearch(), NuFCFitter::GridSearch(), LikelihoodSurfaceBar(), and NuFCDelegateArchiver::UseFitter().

01396 {
01397   mmPars.ContourForNuBars();
01398   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01399 }

Double_t NuSystFitter::FillLikelihoodSurfaceBarGMu ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1405 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNuBarsGMu(), and FillLikelihoodSurfaceGeneric().

01407 {
01408   mmPars.ContourForNuBarsGMu(); //dm2, sn2 and gmu
01409   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01410 }

Double_t NuSystFitter::FillLikelihoodSurfaceCPT ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2522 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), fQuietMode, minimum, NuUtilities::ProgressBar(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

02523 {
02524   Double_t minimum = -1;
02525   const int Nx = surface.GetNbinsX();
02526   const int Ny = surface.GetNbinsY();
02527 
02528   for (Int_t x = 1; x <= Nx; x++)
02529   {
02530       for (Int_t y = 1; y <= Ny; y++)
02531       {
02532           mmPars.Sn2(surface.GetXaxis()->GetBinCenter(x));
02533           mmPars.Dm2(surface.GetYaxis()->GetBinCenter(y));
02534           mmPars.Sn2Bar(surface.GetXaxis()->GetBinCenter(x));
02535           mmPars.Dm2Bar(surface.GetYaxis()->GetBinCenter(y));
02536           Double_t likelihood = (*this)(mmPars.VectorParameters());
02537 
02538           surface.SetBinContent(x, y, likelihood);
02539 
02540           // Work out if this is the minimum
02541           if (minimum < 0 || minimum > likelihood) minimum = likelihood;
02542 
02543           // Now print the progress bar unless in "batch mode"
02544           if(fQuietMode < 2) NuUtilities::ProgressBar((x-1)*Ny+y, Nx*Ny, 1);
02545       }
02546   }
02547 
02548   return minimum;
02549 }

Double_t NuSystFitter::FillLikelihoodSurfaceDecay ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair = 1 
) [virtual]

Definition at line 2492 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDecay(), and FillLikelihoodSurfaceGeneric().

02495 {
02496   mmPars.ContourForDecay(ContourPair);
02497   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02498 }

Double_t NuSystFitter::FillLikelihoodSurfaceDecoherence ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair = 1 
) [virtual]

Definition at line 2502 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDecoherence(), and FillLikelihoodSurfaceGeneric().

02505 {
02506   mmPars.ContourForDecoherence(ContourPair);
02507   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02508 }

Double_t NuSystFitter::FillLikelihoodSurfaceDecoherenceBar ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair = 1 
) [virtual]

Definition at line 2512 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDecoherenceBar(), and FillLikelihoodSurfaceGeneric().

02515 {
02516   mmPars.ContourForDecoherenceBar(ContourPair);
02517   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02518 }

Double_t NuSystFitter::FillLikelihoodSurfaceDm2BarSn2BarGMu3D ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2309 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDm2BarSn2BarGMu3D(), and FillLikelihoodSurfaceGeneric().

02311 {
02312   mmPars.ContourForDm2BarSn2BarGMu3D(); 
02313   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02314 }

Double_t NuSystFitter::FillLikelihoodSurfaceDm2GMuCMu3D ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1440 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDm2GMuCMu3D(), and FillLikelihoodSurfaceGeneric().

01442 {
01443   mmPars.ContourForDm2GMuCMu3D();  
01444   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01445 }

Double_t NuSystFitter::FillLikelihoodSurfaceDm2Sn2GMu3D ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2318 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDm2Sn2GMu3D(), and FillLikelihoodSurfaceGeneric().

02320 {
02321   mmPars.ContourForDm2Sn2GMu3D(); 
02322   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02323 }

Double_t NuSystFitter::FillLikelihoodSurfaceDm2Sn2GMu3DBoth ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2327 of file NuSystFitter.cxx.

References NuMMParameters::ContourForDm2Sn2GMu3DBoth(), and FillLikelihoodSurfaceGeneric().

02329 {
02330   mmPars.ContourForDm2Sn2GMu3DBoth(); 
02331   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02332 }

Double_t NuSystFitter::FillLikelihoodSurfaceGeneric ( TH3D &  surface,
NuMMParameters  mmPars 
) [protected, virtual]

Definition at line 1320 of file NuSystFitter.cxx.

References NuMMParameters::ContourPars3D(), NuMMParameters::FixParameterByIndex(), fQuietMode, Minimise(), minimum, NuUtilities::ProgressBar(), NuMMParameters::SetParameterByIndex(), and NuMMParameters::VectorParameters().

01322 {
01323   const int Nx = surface.GetNbinsX();
01324   const int Ny = surface.GetNbinsY();
01325   const int Nz = surface.GetNbinsZ();
01326 
01327   TAxis* xax = surface.GetXaxis();
01328   TAxis* yax = surface.GetYaxis();
01329   TAxis* zax = surface.GetZaxis();
01330 
01331   //cout << "Before ContourPars3D()" << endl;
01332 
01333   const int idx1 = mmPars.ContourPars3D()[0];
01334   const int idx2 = mmPars.ContourPars3D()[1];
01335   const int idx3 = mmPars.ContourPars3D()[2];
01336 
01337   //cout << "After ContourPars3D(): idx1= " << idx1
01338   //     << ", idx2= " << idx2 << ", idx3= " << idx3 << endl;
01339 
01340   Double_t minimum = -1;
01341 
01342   // Order loops so same-dmsq go sequentially for usual osc contours. This
01343   // takes full advantage of the cacheing in the oscillation functions
01344   for(int z = 1; z <= Nz; z++){
01345     for(int y = 1; y <= Ny; y++){
01346       for(int x = 1; x <= Nx; x++){
01347         // Now print the progress bar unless in "batch mode"
01348         // Need this new mode to avoid having nothing but
01349         // progress bars in FC log files
01350         if(fQuietMode < 2){
01351           NuUtilities::ProgressBar((z-1)*Ny*Nx+y+x-1, Nx*Ny*Nz, 1);
01352         }
01353         //cout << "Before mmPars.SetParameterByIndex()" << endl;
01354         //cout << "(x, y, z): " << x << ", " << y << ", " << z << endl;
01355         //cout << "xax->GetBinCenter(x): " << xax->GetBinCenter(x)
01356         //     << ", yax->GetBinCenter(y): " << yax->GetBinCenter(y)
01357         //     << ", zax->GetBinCenter(z): " << zax->GetBinCenter(z) << endl;
01358         mmPars.SetParameterByIndex(idx3, xax->GetBinCenter(x));
01359         //cout << "Set parameter idx3: " << idx3 << endl;
01360         mmPars.FixParameterByIndex(idx3);
01361         mmPars.SetParameterByIndex(idx2, zax->GetBinCenter(z));
01362         //cout << "Set parameter idx2: " << idx2 << endl;
01363         mmPars.FixParameterByIndex(idx2);
01364         mmPars.SetParameterByIndex(idx1, yax->GetBinCenter(y));
01365         //cout << "Set parameter idx1: " << idx1 << endl;
01366         mmPars.FixParameterByIndex(idx1);
01367         //cout << "In three for loops about to minimize mmPars" << endl;
01368 
01369         mmPars = Minimise(mmPars);
01370 
01371         const Double_t likelihood = (*this)(mmPars.VectorParameters());
01372         surface.SetBinContent(x, y, z, likelihood);
01373 
01374         // Work out if this is the minimum
01375         if (minimum < 0 || minimum > likelihood) minimum = likelihood;
01376       }
01377     }
01378   }
01379   cout << "minimum is: " << minimum << endl;
01380   return minimum;
01381 }

Double_t NuSystFitter::FillLikelihoodSurfaceGeneric ( TH2D &  surface,
NuMMParameters  mmPars,
int  Ix = 1,
int  Iy = 1,
int  Nx = 0,
int  Ny = 0,
TH2D *  mmPars_surfaces[] = 0,
TH2D **  converge_surface = 0,
int  par_index = -1,
int  n_iter = 1,
double *  iter = 0 
) [protected, virtual]

Helper function. Fills a likelihood surface based on the parameters as defined by NuMMParameters::ContourPars()

Definition at line 1229 of file NuSystFitter.cxx.

References NuMMParameters::ContourPars(), fQuietMode, NuMMParameters::GetParameterByIndex(), NuMMParameters::GetParameterName(), NuMMParameters::kNumParameters, Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::SetParameterByIndex().

Referenced by FillLikelihoodSurface(), FillLikelihoodSurfaceBar(), FillLikelihoodSurfaceBarGMu(), FillLikelihoodSurfaceDecay(), FillLikelihoodSurfaceDecoherence(), FillLikelihoodSurfaceDecoherenceBar(), FillLikelihoodSurfaceDm2BarSn2BarGMu3D(), FillLikelihoodSurfaceDm2GMuCMu3D(), FillLikelihoodSurfaceDm2Sn2GMu3D(), FillLikelihoodSurfaceDm2Sn2GMu3DBoth(), FillLikelihoodSurfaceGMuCMuMu(), FillLikelihoodSurfaceGMuGTau(), FillLikelihoodSurfaceGMuGTau3D(), FillLikelihoodSurfaceLED(), FillLikelihoodSurfaceNSI(), FillLikelihoodSurfaceNSI3D(), FillLikelihoodSurfaceNSIBar(), and FillLikelihoodSurfaceSterile().

01234 {
01235 
01236   NuMMParameters surfPars;
01237 
01238   if(!Nx) Nx = surface.GetNbinsX();
01239   if(!Ny) Ny = surface.GetNbinsY();
01240 
01241   TAxis* xax = surface.GetXaxis();
01242   TAxis* yax = surface.GetYaxis();
01243 
01244   const int idx1 = mmPars.ContourPars().first;
01245   const int idx2 = mmPars.ContourPars().second;
01246 
01247   Double_t minimum = -1;
01248 
01249   if(mmPars_surfaces)
01250   {
01251     for(int i = 0; i < NuMMParameters::kNumParameters; i++)
01252     {
01253       const std::string &name = surfPars.GetParameterName(i);
01254       TH2D *clone = (TH2D*)surface.Clone(name.c_str());
01255       clone->SetTitle(name.c_str());
01256       mmPars_surfaces[i] = clone;
01257     }
01258   }
01259 
01260   if(converge_surface) *converge_surface = (TH2D*)surface.Clone("converged");
01261 
01262   // Order loops so same-dmsq go sequentially for usual osc contours. This
01263   // takes full advantage of the cacheing in the oscillation functions
01264   for(int y = Iy; y < Iy+Ny; y++){
01265     for(int x = Ix; x < Ix+Nx; x++){
01266       // Now print the progress bar unless in "batch mode"
01267       // Need this new mode to avoid having nothing but
01268       // progress bars in FC log files
01269       if(fQuietMode < 2){
01270         NuUtilities::ProgressBar((y-Iy)*Nx+(x-Ix+1), Nx*Ny, 1);
01271       }
01272 
01273       double likelihood = 1.0e9;
01274       bool converge = 0;
01275       for(int i = 0; i < n_iter; i++)
01276           {
01277                   NuMMParameters p = mmPars;
01278 
01279                   p.SetParameterByIndex(idx2, xax->GetBinCenter(x));
01280                   p.FixParameterByIndex(idx2);
01281                   p.SetParameterByIndex(idx1, yax->GetBinCenter(y));
01282                   p.FixParameterByIndex(idx1);
01283 
01284           if(par_index >= 0) p.SetParameterByIndex(par_index, iter[i]);
01285 
01286           bool c;
01287                   p = Minimise(p, c);
01288                   double l = (*this)(p.VectorParameters());
01289 
01290           if(l < likelihood)
01291           {
01292             likelihood = l;
01293             converge = c;
01294             surfPars = p;
01295           }
01296           }
01297 
01298       cout << "likelihood: " << likelihood << endl;
01299       surface.SetBinContent(x, y, likelihood);
01300       if(converge_surface) (*converge_surface)->SetBinContent(x,y,converge);
01301 
01302       // Work out if this is the minimum
01303       if (minimum < 0 || minimum > likelihood) minimum = likelihood;
01304 
01305       if(mmPars_surfaces)
01306       {
01307         for(int i = 0; i < NuMMParameters::kNumParameters; i++)
01308         {
01309           mmPars_surfaces[i]->SetBinContent(x, y, surfPars.GetParameterByIndex(i));
01310         }
01311       }
01312     }
01313   }
01314 
01315   return minimum;
01316 }

Double_t NuSystFitter::FillLikelihoodSurfaceGMuCMuMu ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1432 of file NuSystFitter.cxx.

References NuMMParameters::ContourForGMuCMuMu(), and FillLikelihoodSurfaceGeneric().

01434 {                                               
01435   mmPars.ContourForGMuCMuMu();
01436   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01437 }

Double_t NuSystFitter::FillLikelihoodSurfaceGMuGTau ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1414 of file NuSystFitter.cxx.

References NuMMParameters::ContourForGMuGTau(), and FillLikelihoodSurfaceGeneric().

01416 {                                               
01417   mmPars.ContourForGMuGTau();
01418   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01419 }

Double_t NuSystFitter::FillLikelihoodSurfaceGMuGTau3D ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1423 of file NuSystFitter.cxx.

References NuMMParameters::ContourForGMuGTau3D(), and FillLikelihoodSurfaceGeneric().

01425 {
01426   mmPars.ContourForGMuGTau3D(); // 3rd parameter is cmumu
01427   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01428 }

Double_t NuSystFitter::FillLikelihoodSurfaceLED ( TH2D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2565 of file NuSystFitter.cxx.

References NuMMParameters::ContourForLED(), and FillLikelihoodSurfaceGeneric().

02565                                                                                    {
02566   mmPars.ContourForLED();
02567   return FillLikelihoodSurfaceGeneric(surface, mmPars);    
02568 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSI ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair 
) [virtual]

Definition at line 2336 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNSI(), and FillLikelihoodSurfaceGeneric().

Referenced by NuFCFitterNSINu::Archive(), NuFCFitterNSINu::GridSearch(), and GridSearchNSI().

02339 {
02340   mmPars.ContourForNSI(ContourPair);
02341   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02342 
02343 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSI3D ( TH3D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2356 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNSI(), and FillLikelihoodSurfaceGeneric().

Referenced by NuFCFitterNSI::Archive(), and NuFCFitterNSI::GridSearch().

02358 {
02359   mmPars.ContourForNSI();
02360   //  cout << "Set contours for NSI: ContourPars3D()[0]= " << mmPars.ContourPars3D()[0]
02361   //     << ", ContourPars3D()[1]= " << mmPars.ContourPars3D()[1]
02362   //     << ", ContourPars3D()[2]= " << mmPars.ContourPars3D()[2]
02363   //     << endl;
02364   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02365 
02366 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSIBar ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair 
) [virtual]

Definition at line 2346 of file NuSystFitter.cxx.

References NuMMParameters::ContourForNSIBar(), and FillLikelihoodSurfaceGeneric().

Referenced by NuFCFitterNSINubar::Archive(), and NuFCFitterNSINubar::GridSearch().

02349 {
02350   mmPars.ContourForNSIBar(ContourPair);
02351   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02352 
02353 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSIDm2Eps ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  sinvals 
) [virtual]

Definition at line 2371 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Epsilon(), NuMMParameters::FixDm2(), NuMMParameters::FixEpsilon(), fQuietMode, Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

02374 {
02375  
02376   Int_t Neps = surface.GetNbinsX();
02377   Int_t Ndm2 = surface.GetNbinsY();
02378 
02379   Double_t minimum = 1e6;
02380 
02381   for (int i = 1; i <= Neps; i++){
02382     for (int j = 1; j <= Ndm2; j++){
02383 
02384       if(fQuietMode < 2){
02385         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Neps, 1);
02386       }
02387 
02388       Double_t eps = surface.GetXaxis()->GetBinCenter(i);
02389       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
02390       mmPars.Dm2(dm2);
02391       mmPars.Epsilon(eps);
02392       mmPars.FixDm2();
02393       mmPars.FixEpsilon();
02394 
02395       NuMMParameters bf = this->Minimise(mmPars);
02396       Double_t likelihood = Likelihood(bf);
02397       Double_t sn2Min = bf.Sn2();
02398 
02399       surface.SetBinContent(i, j, likelihood);
02400       sinvals.SetBinContent(i, j, sn2Min);
02401       if(likelihood < minimum) minimum = likelihood;      
02402     }
02403   }
02404 
02405   return minimum;
02406 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSIDm2Sn2 ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  epsvals 
) [virtual]

Definition at line 2451 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Epsilon(), NuMMParameters::FixDm2(), NuMMParameters::FixSn2(), fQuietMode, Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

02454 {
02455  
02456   Int_t Nsn2 = surface.GetNbinsX();
02457   Int_t Ndm2 = surface.GetNbinsY();
02458 
02459   Double_t minimum = 1e6;
02460 
02461   for (int i = 1; i <= Nsn2; i++){
02462     for (int j = 1; j <= Ndm2; j++){
02463 
02464       if(fQuietMode < 2){
02465         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Nsn2, 1);
02466       }
02467 
02468       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
02469       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
02470       mmPars.Dm2(dm2);
02471       mmPars.Sn2(sn2);
02472       mmPars.FixDm2();
02473       mmPars.FixSn2();
02474 
02475       NuMMParameters bf = this->Minimise(mmPars);
02476       Double_t likelihood = Likelihood(bf);
02477       Double_t epsMin = bf.Epsilon();
02478       
02479       surface.SetBinContent(i, j, likelihood);
02480       epsvals.SetBinContent(i, j, epsMin);
02481 
02482       if(likelihood < minimum) minimum = likelihood;      
02483     }
02484   }
02485 
02486   return minimum;
02487 }

Double_t NuSystFitter::FillLikelihoodSurfaceNSISn2Eps ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  massvals 
) [virtual]

Definition at line 2410 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Epsilon(), NuMMParameters::FixEpsilon(), NuMMParameters::FixSn2(), fQuietMode, Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

02413 {
02414  
02415   Int_t Neps = surface.GetNbinsX();
02416   Int_t Nsn2 = surface.GetNbinsY();
02417 
02418   Double_t minimum = 1e6;
02419 
02420   for (int i = 1; i <= Neps; i++){
02421     for (int j = 1; j <= Nsn2; j++){
02422 
02423       if(fQuietMode < 2){
02424         NuUtilities::ProgressBar((i-1)*Nsn2+j, Nsn2*Neps, 1);
02425       }
02426 
02427       Double_t eps = surface.GetXaxis()->GetBinCenter(i);
02428       Double_t sn2 = surface.GetYaxis()->GetBinCenter(j);
02429       mmPars.Sn2(sn2);
02430       mmPars.Epsilon(eps);
02431       mmPars.FixSn2();
02432       mmPars.FixEpsilon();
02433 
02434       NuMMParameters bf = this->Minimise(mmPars);
02435       Double_t likelihood = Likelihood(bf);
02436       Double_t dm2Min = bf.Dm2();
02437 
02438       surface.SetBinContent(i, j, likelihood);
02439       massvals.SetBinContent(i, j, dm2Min);
02440 
02441       if(likelihood < minimum) minimum = likelihood;      
02442     }
02443   }
02444 
02445   return minimum;
02446 }

Double_t NuSystFitter::FillLikelihoodSurfaceSterile ( TH2D &  surface,
NuMMParameters  mmPars,
Int_t  ContourPair,
int  Ix = 1,
int  Iy = 1,
int  Nx = 0,
int  Ny = 0,
TH2D *  mmPars_surfaces[] = 0,
TH2D **  converge_surface = 0,
int  par_index = -1,
int  n_iter = 1,
double *  iter = 0 
) [virtual]

Definition at line 2552 of file NuSystFitter.cxx.

References NuMMParameters::ContourForSterile(), and FillLikelihoodSurfaceGeneric().

02557 {
02558   mmPars.ContourForSterile(ContourPair);
02559   return FillLikelihoodSurfaceGeneric(surface, mmPars, Ix, Iy, Nx, Ny, mmPars_surfaces, converge_surface, par_index, n_iter, iter);
02560 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransCMuMuBoth ( TH1D &  surface,
TH1D &  gmuvals,
TH1D &  dm2vals,
TH1D &  sn2vals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1943 of file NuSystFitter.cxx.

References NuMMParameters::CMuMu(), NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2Bar(), NuMMParameters::FixCMuMu(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseSn2(), NuMMParameters::ReleaseSn2Bar(), and NuMMParameters::Sn2Bar().

01948 {
01949 
01950   Int_t Ncmumu = surface.GetNbinsX();
01951   Double_t minimum = 1e6;
01952   for (int i = 1; i <= Ncmumu; i++){
01953 
01954     //  if(fQuietMode < 2){
01955     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ncmumu, 1);
01956     //  }
01957 
01958       Double_t cmumu = surface.GetXaxis()->GetBinCenter(i);
01959       mmPars.CMuMu(cmumu);
01960       mmPars.FixCMuMu();
01961 
01962       mmPars.ReleaseDm2Bar();
01963       mmPars.ReleaseSn2Bar();
01964       mmPars.ReleaseDm2();
01965       mmPars.ReleaseSn2();
01966       mmPars.ConstrainCPT();
01967 
01968       NuMMParameters bf = this->Minimise(mmPars);
01969       Double_t likelihood = Likelihood(bf);
01970       Double_t gmuMin = bf.GMu();
01971       Double_t sn2Min = bf.Sn2Bar();
01972       Double_t dm2Min = bf.Dm2Bar();
01973 
01974       surface.SetBinContent(i, likelihood);
01975       gmuvals.SetBinContent(i, gmuMin);
01976       dm2vals.SetBinContent(i, dm2Min);
01977       sn2vals.SetBinContent(i, sn2Min);
01978       if(likelihood < minimum) minimum = likelihood;
01979   }
01980   return minimum;
01981 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransCTauTauBoth ( TH1D &  surface,
TH1D &  gmuvals,
TH1D &  dm2vals,
TH1D &  sn2vals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1984 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainCPT(), NuMMParameters::CTauTau(), NuMMParameters::Dm2Bar(), NuMMParameters::FixCTauTau(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseSn2(), NuMMParameters::ReleaseSn2Bar(), and NuMMParameters::Sn2Bar().

01989 {                                                          
01990 
01991   Int_t Nctautau = surface.GetNbinsX();
01992   Double_t minimum = 1e6;
01993   for (int i = 1; i <= Nctautau; i++){
01994 
01995     //  if(fQuietMode < 2){
01996     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Nctautau, 1);
01997     //  }
01998 
01999       Double_t ctautau = surface.GetXaxis()->GetBinCenter(i);
02000       mmPars.CTauTau(ctautau);
02001       mmPars.FixCTauTau();
02002 
02003       mmPars.ReleaseDm2Bar();
02004       mmPars.ReleaseSn2Bar();
02005       mmPars.ReleaseDm2();
02006       mmPars.ReleaseSn2();
02007       mmPars.ConstrainCPT();
02008 
02009       NuMMParameters bf = this->Minimise(mmPars);
02010       Double_t likelihood = Likelihood(bf);
02011       Double_t gmuMin = bf.GMu();
02012       Double_t sn2Min = bf.Sn2Bar();
02013       Double_t dm2Min = bf.Dm2Bar();
02014 
02015       surface.SetBinContent(i, likelihood);
02016       gmuvals.SetBinContent(i, gmuMin);
02017       dm2vals.SetBinContent(i, dm2Min);
02018       sn2vals.SetBinContent(i, sn2Min);
02019       if(likelihood < minimum) minimum = likelihood;
02020   }
02021   return minimum;
02022 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2BarGMu ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  sinvals 
) [virtual]

Definition at line 1449 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixGMu(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2Bar().

01452 {
01453 
01454   Int_t Ngmu = surface.GetNbinsX();
01455   Int_t Ndm2 = surface.GetNbinsY();
01456 
01457   Double_t minimum = 1e6;
01458 
01459   for (int i = 1; i <= Ngmu; i++){
01460     for (int j = 1; j <= Ndm2; j++){
01461 
01462       if(fQuietMode < 2){
01463         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01464       }
01465 
01466       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01467       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01468       mmPars.Dm2Bar(dm2);
01469       mmPars.GMu(gmu);
01470       mmPars.FixDm2Bar();
01471       mmPars.FixGMu();
01472 
01473       NuMMParameters bf = this->Minimise(mmPars);
01474       Double_t likelihood = Likelihood(bf);
01475       Double_t sn2Min = bf.Sn2Bar();
01476 
01477       surface.SetBinContent(i, j, likelihood);
01478       sinvals.SetBinContent(i, j, sn2Min);
01479       if(likelihood < minimum) minimum = likelihood;
01480     }
01481   }
01482 
01483   return minimum;
01484 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2BarSn2Bar ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  gmuvals 
) [virtual]

Definition at line 1686 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixSn2Bar(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2Bar().

01689 {
01690 
01691   Int_t Nsn2 = surface.GetNbinsX();
01692   Int_t Ndm2 = surface.GetNbinsY();
01693 
01694   Double_t minimum = 1e6;
01695 
01696   for (int i = 1; i <= Nsn2; i++){
01697     for (int j = 1; j <= Ndm2; j++){
01698 
01699       if(fQuietMode < 2){
01700         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Nsn2, 1);
01701       }
01702 
01703       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
01704       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01705       mmPars.Dm2Bar(dm2);
01706       mmPars.Sn2Bar(sn2);
01707       mmPars.FixDm2Bar();
01708       mmPars.FixSn2Bar();
01709 
01710       NuMMParameters bf = this->Minimise(mmPars);
01711       Double_t likelihood = Likelihood(bf);
01712       Double_t BMin = bf.GMu();
01713 
01714       surface.SetBinContent(i, j, likelihood);
01715       bvals.SetBinContent(i, j, BMin);
01716 
01717       if(likelihood < minimum) minimum = likelihood;
01718     }
01719   }
01720 
01721   return minimum;
01722 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2Both ( TH1D &  surface,
TH1D &  sn2vals,
TH1D &  gmuvals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1855 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseGMu(), NuMMParameters::ReleaseSn2(), NuMMParameters::ReleaseSn2Bar(), and NuMMParameters::Sn2Bar().

01859 {
01860 
01861   Int_t Ndm2 = surface.GetNbinsX();
01862 
01863   Double_t minimum = 1e6;
01864 
01865   for (int i = 1; i <= Ndm2; i++){
01866 
01867     //  if(fQuietMode < 2){
01868     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01869     //  }
01870 
01871       Double_t dm2 = surface.GetXaxis()->GetBinCenter(i);
01872       mmPars.Dm2(dm2);
01873       mmPars.Dm2Bar(dm2);
01874       mmPars.FixDm2();   
01875       mmPars.FixDm2Bar();   
01876 
01877       mmPars.ReleaseSn2();
01878       mmPars.ReleaseSn2Bar();
01879       mmPars.ReleaseGMu();
01880       mmPars.ConstrainCPT();
01881 
01882       NuMMParameters bf = this->Minimise(mmPars);
01883       Double_t likelihood = Likelihood(bf);
01884       Double_t sn2Min = bf.Sn2Bar();
01885       Double_t gmuMin = bf.GMu();
01886       cout<<"i="<<i<<", dm2="<<dm2<<", likelihood="<<likelihood<<endl;
01887 
01888 
01889       surface.SetBinContent(i, likelihood);
01890       sn2vals.SetBinContent(i,sn2Min);
01891       gmuvals.SetBinContent(i,gmuMin);
01892       if(likelihood < minimum) minimum = likelihood;
01893   }
01894 
01895   return minimum;
01896 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2GMu ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  sinvals 
) [virtual]

Definition at line 1648 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::FixDm2(), NuMMParameters::FixGMu(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

01651 {
01652 
01653   Int_t Ngmu = surface.GetNbinsX();
01654   Int_t Ndm2 = surface.GetNbinsY();
01655 
01656   Double_t minimum = 1e6;
01657 
01658   for (int i = 1; i <= Ngmu; i++){
01659     for (int j = 1; j <= Ndm2; j++){
01660 
01661       if(fQuietMode < 2){
01662         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01663       }
01664 
01665       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01666       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01667       mmPars.Dm2(dm2);
01668       mmPars.GMu(gmu);
01669       mmPars.FixDm2();
01670       mmPars.FixGMu();
01671 
01672       NuMMParameters bf = this->Minimise(mmPars);
01673       Double_t likelihood = Likelihood(bf);
01674       Double_t sn2Min = bf.Sn2();
01675 
01676       surface.SetBinContent(i, j, likelihood);
01677       sinvals.SetBinContent(i, j, sn2Min);
01678       if(likelihood < minimum) minimum = likelihood;
01679     }
01680   }
01681 
01682   return minimum;
01683 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2GMuBoth ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  sinvals 
) [virtual]

Definition at line 1486 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixGMu(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2Bar().

01489 { 
01490   
01491   Int_t Ngmu = surface.GetNbinsX();
01492   Int_t Ndm2 = surface.GetNbinsY();
01493 
01494   Double_t minimum = 1e6;
01495                                             
01496   for (int i = 1; i <= Ngmu; i++){
01497     for (int j = 1; j <= Ndm2; j++){
01498   
01499       if(fQuietMode < 2){
01500         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01501       }
01502 
01503       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01504       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01505       mmPars.Dm2Bar(dm2);
01506       mmPars.Dm2(dm2);
01507       mmPars.GMu(gmu);
01508       mmPars.FixDm2Bar();
01509       mmPars.FixDm2();
01510       mmPars.ConstrainCPT();
01511       mmPars.FixGMu();
01512 
01513       NuMMParameters bf = this->Minimise(mmPars);
01514       Double_t likelihood = Likelihood(bf);
01515       Double_t sn2Min = bf.Sn2Bar();
01516 
01517       surface.SetBinContent(i, j, likelihood);
01518       sinvals.SetBinContent(i, j, sn2Min);
01519       if(likelihood < minimum) minimum = likelihood;
01520     }
01521   }
01522 
01523   return minimum;
01524 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2NQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2202 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::FixDm2(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseGMu(), and NuMMParameters::ReleaseSn2().

02204 {
02205 
02206   Int_t Ndm2 = surface.GetNbinsX();
02207 
02208   Double_t minimum = 1e6;
02209 
02210   for (int i = 1; i <= Ndm2; i++){
02211 
02212     //  if(fQuietMode < 2){
02213     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02214     //  }
02215 
02216       Double_t dm2 = surface.GetXaxis()->GetBinCenter(i);
02217       mmPars.Dm2(dm2);
02218       mmPars.FixDm2();   
02219 
02220       mmPars.ReleaseSn2();
02221       mmPars.ReleaseGMu();
02222 
02223       NuMMParameters bf = this->Minimise(mmPars);
02224       Double_t likelihood = Likelihood(bf);
02225 
02226       surface.SetBinContent(i, likelihood);
02227       if(likelihood < minimum) minimum = likelihood;
02228   }
02229 
02230   return minimum;
02231 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2PQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2101 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2Bar(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseGMu(), and NuMMParameters::ReleaseSn2Bar().

02103 {
02104 
02105   Int_t Ndm2 = surface.GetNbinsX();
02106 
02107   Double_t minimum = 1e6;
02108 
02109   for (int i = 1; i <= Ndm2; i++){
02110 
02111     //  if(fQuietMode < 2){
02112     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02113     //  }
02114 
02115       Double_t dm2 = surface.GetXaxis()->GetBinCenter(i);
02116       mmPars.Dm2Bar(dm2);
02117       mmPars.FixDm2Bar();   
02118 
02119       mmPars.ReleaseSn2Bar();
02120       mmPars.ReleaseGMu();
02121 
02122       NuMMParameters bf = this->Minimise(mmPars);
02123       Double_t likelihood = Likelihood(bf);
02124 
02125       surface.SetBinContent(i, likelihood);
02126       if(likelihood < minimum) minimum = likelihood;
02127   }
02128 
02129   return minimum;
02130 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2Sn2 ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  gmuvals 
) [virtual]

Definition at line 1770 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::FixDm2(), NuMMParameters::FixSn2(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

01773 {
01774 
01775   Int_t Nsn2 = surface.GetNbinsX();
01776   Int_t Ndm2 = surface.GetNbinsY();
01777 
01778   Double_t minimum = 1e6;
01779 
01780   for (int i = 1; i <= Nsn2; i++){
01781     for (int j = 1; j <= Ndm2; j++){
01782 
01783       if(fQuietMode < 2){
01784         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Nsn2, 1);
01785       }
01786 
01787       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
01788       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01789       mmPars.Dm2(dm2);
01790       mmPars.Sn2(sn2);
01791       mmPars.FixDm2();
01792       mmPars.FixSn2();
01793 
01794       NuMMParameters bf = this->Minimise(mmPars);
01795       Double_t likelihood = Likelihood(bf);
01796       Double_t BMin = bf.GMu();
01797 
01798       surface.SetBinContent(i, j, likelihood);
01799       bvals.SetBinContent(i, j, BMin);
01800 
01801       if(likelihood < minimum) minimum = likelihood;
01802     }
01803   }
01804 
01805   return minimum;
01806 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransDm2Sn2Both ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  gmuvals 
) [virtual]

Definition at line 1725 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::FixDm2(), NuMMParameters::FixDm2Bar(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

01728 {                                                            
01729 
01730   Int_t Nsn2 = surface.GetNbinsX();
01731   Int_t Ndm2 = surface.GetNbinsY();
01732 
01733   Double_t minimum = 1e6;
01734 
01735   for (int i = 1; i <= Nsn2; i++){
01736     for (int j = 1; j <= Ndm2; j++){
01737 
01738       if(fQuietMode < 2){
01739         NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Nsn2, 1);
01740       }
01741 
01742       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
01743       Double_t dm2 = surface.GetYaxis()->GetBinCenter(j);
01744       mmPars.Dm2Bar(dm2);
01745       mmPars.Sn2Bar(sn2);
01746       mmPars.FixDm2Bar();
01747       mmPars.FixSn2Bar();
01748       mmPars.Dm2(dm2);
01749       mmPars.Sn2(sn2);
01750       mmPars.FixDm2();
01751       mmPars.FixSn2();
01752 
01753       NuMMParameters bf = this->Minimise(mmPars);
01754       Double_t likelihood = Likelihood(bf);
01755       Double_t BMin = bf.GMu();
01756 
01757       surface.SetBinContent(i, j, likelihood);
01758       bvals.SetBinContent(i, j, BMin);
01759 
01760       if(likelihood < minimum) minimum = likelihood;
01761     }
01762   }
01763 
01764   return minimum;
01765 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransGMuBoth ( TH1D &  surface,
TH1D &  cmuvals,
TH1D &  dm2vals,
TH1D &  sn2vals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1810 of file NuSystFitter.cxx.

References NuMMParameters::CMuMu(), NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2Bar(), NuMMParameters::FixGMu(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseSn2(), NuMMParameters::ReleaseSn2Bar(), and NuMMParameters::Sn2Bar().

01815 {     
01816       
01817   Int_t Ngmu = surface.GetNbinsX();
01818 
01819   Double_t minimum = 1e6; 
01820       
01821   for (int i = 1; i <= Ngmu; i++){
01822       
01823     //  if(fQuietMode < 2){
01824     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01825     //  }
01826     
01827       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01828       mmPars.GMu(gmu);
01829       mmPars.FixGMu();
01830 
01831       mmPars.ReleaseDm2Bar();
01832       mmPars.ReleaseSn2Bar();
01833       mmPars.ReleaseDm2();
01834       mmPars.ReleaseSn2();
01835       mmPars.ConstrainCPT();
01836 
01837       NuMMParameters bf = this->Minimise(mmPars);
01838       Double_t likelihood = Likelihood(bf);
01839       Double_t cmuMin = bf.CMuMu();
01840       Double_t sn2Min = bf.Sn2Bar();
01841       Double_t dm2Min = bf.Dm2Bar();
01842 
01843       surface.SetBinContent(i, likelihood);
01844       cmuvals.SetBinContent(i, cmuMin);
01845       dm2vals.SetBinContent(i, dm2Min);
01846       sn2vals.SetBinContent(i, sn2Min);
01847       if(likelihood < minimum) minimum = likelihood;
01848   }
01849 
01850   return minimum;
01851 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransGMuNQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2169 of file NuSystFitter.cxx.

References NuMMParameters::FixGMu(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), and NuMMParameters::ReleaseSn2().

02171 {     
02172       
02173   Int_t Ngmu = surface.GetNbinsX();
02174 
02175   Double_t minimum = 1e6; 
02176       
02177   for (int i = 1; i <= Ngmu; i++){
02178       
02179     //  if(fQuietMode < 2){
02180     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02181     //  }
02182     
02183       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
02184       mmPars.GMu(gmu);
02185       mmPars.FixGMu();
02186 
02187       mmPars.ReleaseDm2();
02188       mmPars.ReleaseSn2();
02189 
02190       NuMMParameters bf = this->Minimise(mmPars);
02191       Double_t likelihood = Likelihood(bf);
02192 
02193       surface.SetBinContent(i, likelihood);
02194       if(likelihood < minimum) minimum = likelihood;
02195   }
02196 
02197   return minimum;
02198 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransGMuPQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2068 of file NuSystFitter.cxx.

References NuMMParameters::FixGMu(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2Bar(), and NuMMParameters::ReleaseSn2Bar().

02070 {     
02071       
02072   Int_t Ngmu = surface.GetNbinsX();
02073 
02074   Double_t minimum = 1e6; 
02075       
02076   for (int i = 1; i <= Ngmu; i++){
02077       
02078     //  if(fQuietMode < 2){
02079     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02080     //  }
02081     
02082       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
02083       mmPars.GMu(gmu);
02084       mmPars.FixGMu();
02085 
02086       mmPars.ReleaseDm2Bar();
02087       mmPars.ReleaseSn2Bar();
02088 
02089       NuMMParameters bf = this->Minimise(mmPars);
02090       Double_t likelihood = Likelihood(bf);
02091 
02092       surface.SetBinContent(i, likelihood);
02093       if(likelihood < minimum) minimum = likelihood;
02094   }
02095 
02096   return minimum;
02097 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransGTauBoth ( TH1D &  surface,
TH1D &  cmuvals,
TH1D &  dm2vals,
TH1D &  sn2vals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2024 of file NuSystFitter.cxx.

References NuMMParameters::CMuMu(), NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2Bar(), NuMMParameters::FixGTau(), NuMMParameters::GTau(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseSn2(), NuMMParameters::ReleaseSn2Bar(), and NuMMParameters::Sn2Bar().

02029 {                                                          
02030 
02031   Int_t Ngtau = surface.GetNbinsX();
02032   Double_t minimum = 1e6;
02033   for (int i = 1; i <= Ngtau; i++){
02034 
02035     //  if(fQuietMode < 2){
02036     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngtau, 1);
02037     //  }
02038 
02039       Double_t gtau = surface.GetXaxis()->GetBinCenter(i);
02040       mmPars.GTau(gtau);
02041       mmPars.FixGTau();
02042 
02043       mmPars.ReleaseDm2Bar();
02044       mmPars.ReleaseSn2Bar();
02045       mmPars.ReleaseDm2();
02046       mmPars.ReleaseSn2();
02047       mmPars.ConstrainCPT();
02048 
02049       NuMMParameters bf = this->Minimise(mmPars);
02050       Double_t likelihood = Likelihood(bf);
02051       Double_t cmuMin = bf.CMuMu();
02052       Double_t sn2Min = bf.Sn2Bar();
02053       Double_t dm2Min = bf.Dm2Bar();
02054 
02055       surface.SetBinContent(i, likelihood);
02056       cmuvals.SetBinContent(i, cmuMin);
02057       dm2vals.SetBinContent(i, dm2Min);
02058       sn2vals.SetBinContent(i, sn2Min);
02059       if(likelihood < minimum) minimum = likelihood;
02060   }
02061   return minimum;
02062 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2BarGMu ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  massvals 
) [virtual]

Definition at line 1528 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixGMu(), NuMMParameters::FixSn2Bar(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2Bar().

01531 {
01532 
01533   Int_t Ngmu = surface.GetNbinsX();
01534   Int_t Nsn2 = surface.GetNbinsY();
01535 
01536   Double_t minimum = 1e6;
01537 
01538   for (int i = 1; i <= Ngmu; i++){
01539     for (int j = 1; j <= Nsn2; j++){
01540 
01541       if(fQuietMode < 2){
01542         NuUtilities::ProgressBar((i-1)*Nsn2+j, Nsn2*Ngmu, 1);
01543       }
01544 
01545       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01546       Double_t sn2 = surface.GetYaxis()->GetBinCenter(j);
01547       mmPars.Sn2Bar(sn2);
01548       mmPars.GMu(gmu);
01549       mmPars.FixSn2Bar();
01550       mmPars.FixGMu();
01551 
01552       NuMMParameters bf = this->Minimise(mmPars);
01553       Double_t likelihood = Likelihood(bf);
01554       Double_t dm2Min = bf.Dm2Bar();
01555 
01556       surface.SetBinContent(i, j, likelihood);
01557       massvals.SetBinContent(i, j, dm2Min);
01558       if(likelihood < minimum) minimum = likelihood;
01559     }
01560   }
01561 
01562   return minimum;
01563 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2Both ( TH1D &  surface,
TH1D &  dm2vals,
TH1D &  gmuvals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 1900 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2Bar(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseGMu(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

01904 {
01905 
01906   Int_t Nsn2 = surface.GetNbinsX();
01907 
01908   Double_t minimum = 1e6;
01909 
01910   for (int i = 1; i <= Nsn2; i++){                
01911                                                   
01912     //  if(fQuietMode < 2){
01913     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
01914     //  }
01915 
01916       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
01917       mmPars.Sn2(sn2);
01918       mmPars.Sn2Bar(sn2);
01919       mmPars.FixSn2();   
01920       mmPars.FixSn2Bar();                        
01921                                                   
01922       mmPars.ReleaseDm2();
01923       mmPars.ReleaseDm2Bar();
01924       mmPars.ReleaseGMu();
01925       mmPars.ConstrainCPT();
01926 
01927       NuMMParameters bf = this->Minimise(mmPars);
01928       Double_t likelihood = Likelihood(bf);
01929       Double_t dm2Min = bf.Dm2Bar();
01930       Double_t gmuMin = bf.GMu();                
01931                                                   
01932       surface.SetBinContent(i, likelihood);
01933       dm2vals.SetBinContent(i, dm2Min);
01934       gmuvals.SetBinContent(i, gmuMin);
01935       if(likelihood < minimum) minimum = likelihood;
01936   }
01937 
01938   return minimum;
01939 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2GMu ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  massvals 
) [virtual]

Definition at line 1608 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::FixGMu(), NuMMParameters::FixSn2(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), and NuMMParameters::Sn2().

01611 {
01612 
01613   Int_t Ngmu = surface.GetNbinsX();
01614   Int_t Nsn2 = surface.GetNbinsY();
01615 
01616   Double_t minimum = 1e6;
01617 
01618   for (int i = 1; i <= Ngmu; i++){
01619     for (int j = 1; j <= Nsn2; j++){
01620 
01621       if(fQuietMode < 2){
01622         NuUtilities::ProgressBar((i-1)*Nsn2+j, Nsn2*Ngmu, 1);
01623       }
01624 
01625       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01626       Double_t sn2 = surface.GetYaxis()->GetBinCenter(j);
01627       mmPars.Sn2(sn2);
01628       mmPars.GMu(gmu);
01629       mmPars.FixSn2();
01630       mmPars.FixGMu();
01631 
01632       NuMMParameters bf = this->Minimise(mmPars);
01633       Double_t likelihood = Likelihood(bf);
01634       Double_t dm2Min = bf.Dm2();
01635 
01636       surface.SetBinContent(i, j, likelihood);
01637       massvals.SetBinContent(i, j, dm2Min);
01638       if(likelihood < minimum) minimum = likelihood;
01639     }
01640   }
01641 
01642   return minimum;
01643 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2GMuBoth ( TH2D &  surface,
NuMMParameters  mmPars,
TH2D &  massvals 
) [virtual]

Definition at line 1565 of file NuSystFitter.cxx.

References NuMMParameters::ConstrainCPT(), NuMMParameters::Dm2Bar(), NuMMParameters::FixGMu(), NuMMParameters::FixSn2(), NuMMParameters::FixSn2Bar(), fQuietMode, NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, NuUtilities::ProgressBar(), NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

01568 {
01569 
01570   Int_t Ngmu = surface.GetNbinsX();
01571   Int_t Nsn2 = surface.GetNbinsY();
01572 
01573   Double_t minimum = 1e6;
01574 
01575   for (int i = 1; i <= Ngmu; i++){
01576     for (int j = 1; j <= Nsn2; j++){
01577 
01578       if(fQuietMode < 2){
01579         NuUtilities::ProgressBar((i-1)*Nsn2+j, Nsn2*Ngmu, 1);
01580       }
01581       
01582       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
01583       Double_t sn2 = surface.GetYaxis()->GetBinCenter(j);
01584       mmPars.Sn2Bar(sn2);
01585       mmPars.Sn2(sn2);
01586       mmPars.GMu(gmu);
01587       mmPars.FixSn2Bar(); 
01588       mmPars.FixSn2();
01589       mmPars.ConstrainCPT();
01590       mmPars.FixGMu();
01591 
01592       NuMMParameters bf = this->Minimise(mmPars);
01593       Double_t likelihood = Likelihood(bf);
01594       Double_t dm2Min = bf.Dm2Bar();
01595     
01596       surface.SetBinContent(i, j, likelihood);
01597       massvals.SetBinContent(i, j, dm2Min);
01598       if(likelihood < minimum) minimum = likelihood;
01599     }
01600   }
01601 
01602   return minimum;
01603 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2NQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2235 of file NuSystFitter.cxx.

References NuMMParameters::FixSn2(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2(), NuMMParameters::ReleaseGMu(), and NuMMParameters::Sn2().

02237 {
02238 
02239   Int_t Nsn2 = surface.GetNbinsX();
02240 
02241   Double_t minimum = 1e6;
02242 
02243   for (int i = 1; i <= Nsn2; i++){                
02244                                                   
02245     //  if(fQuietMode < 2){
02246     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02247     //  }
02248 
02249       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
02250       mmPars.Sn2(sn2);
02251       mmPars.FixSn2();                        
02252                                                   
02253       mmPars.ReleaseDm2();
02254       mmPars.ReleaseGMu();
02255 
02256       NuMMParameters bf = this->Minimise(mmPars);
02257       Double_t likelihood = Likelihood(bf);
02258                                                   
02259       surface.SetBinContent(i, likelihood);
02260       if(likelihood < minimum) minimum = likelihood;
02261   }
02262 
02263   return minimum;
02264 }

Double_t NuSystFitter::FillLikelihoodSurfaceTransSn2PQ ( TH1D &  surface,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2134 of file NuSystFitter.cxx.

References NuMMParameters::FixSn2Bar(), Likelihood(), Minimise(), minimum, NuMMParameters::ReleaseDm2Bar(), NuMMParameters::ReleaseGMu(), and NuMMParameters::Sn2Bar().

02136 {
02137 
02138   Int_t Nsn2 = surface.GetNbinsX();
02139 
02140   Double_t minimum = 1e6;
02141 
02142   for (int i = 1; i <= Nsn2; i++){                
02143                                                   
02144     //  if(fQuietMode < 2){
02145     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02146     //  }
02147 
02148       Double_t sn2 = surface.GetXaxis()->GetBinCenter(i);
02149       mmPars.Sn2Bar(sn2);
02150       mmPars.FixSn2Bar();                        
02151                                                   
02152       mmPars.ReleaseDm2Bar();
02153       mmPars.ReleaseGMu();
02154 
02155       NuMMParameters bf = this->Minimise(mmPars);
02156       Double_t likelihood = Likelihood(bf);
02157                                                   
02158       surface.SetBinContent(i, likelihood);
02159       if(likelihood < minimum) minimum = likelihood;
02160   }
02161 
02162   return minimum;
02163 }

NuMMParameters NuSystFitter::FineGridSearch ( NuMMParameters  mmStart  )  [protected]

Definition at line 2639 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), MuELoss::e, GridSearch(), Msg::kInfo, NuMMParameters::LowerDm2BarLimit(), NuMMParameters::LowerSn2BarLimit(), MSG, and NuMMParameters::Sn2Bar().

Referenced by FCMinimise().

02640 {
02641   // Calculate the location of the fine grid search
02642   Double_t gridL, gridR, gridT, gridB;
02643   gridL = gridR = gridT = gridB = -1.0;
02644 
02645   // Go +- 0.2 in sin space
02646   gridL = mmStart.Sn2Bar() - 0.2;
02647   gridR = gridL + 0.4;
02648   // Now limit it to the physical realm
02649   if (gridL <= mmStart.LowerSn2BarLimit()) gridL = mmStart.LowerSn2BarLimit();
02650   if (gridR > 1.0) gridR = 1.0;
02651 
02652   // And go +-2e-3 in mass space
02653   gridB = mmStart.Dm2Bar() - 2e-3;
02654   gridT = gridB + 4e-3;
02655   // We only really need to worry about lower limit here
02656   if (gridB <= mmStart.LowerDm2BarLimit()) gridB = mmStart.LowerSn2BarLimit();
02657 
02658 
02659   MSG("NuSystFitter", Msg::kInfo)
02660     << "Doing fine grid search on sin = [ "
02661     << gridL << ", " << gridR << " ]\n"
02662     << "                              dm  = [ "
02663     << gridB << ", " << gridT << " ]\n";
02664 
02665   // Ensure that we haven't done anything stupid like set low above high
02666   if (gridL >= gridR || gridT <= gridB) {
02667 //      cout << "Error: Negative range in fine grid search" << endl;
02668       throw runtime_error("Negative range in fine grid search");
02669   }
02670 
02671   // Make the surface
02672 
02673   // The width and height of the grid
02674   const Int_t gridWidth = 19;
02675   const Int_t gridHeight = 19;
02676 
02677   // Generate the surface
02678   TH2D likesurf("likesurf", "likesurf", gridWidth, gridL, gridR,
02679     gridHeight, gridB, gridT);
02680 
02681   // And now do a grid search on this
02682   NuMMParameters mmGrid = GridSearch(likesurf, mmStart);
02683 
02684   // Output the minimum
02685   MSG("NusystFitter", Msg::kInfo) << "Fine grid search found sn2bar: "
02686     << mmGrid.Sn2Bar() << ", dm2bar: " << mmGrid.Dm2Bar() << '\n';
02687 
02688   return mmGrid;
02689 }

NuMMParameters NuSystFitter::FineGridSearchNSI ( NuMMParameters  mmStart  )  [protected]

Definition at line 2959 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), MuELoss::e, NuMMParameters::Epsilon(), GridSearchNSI(), Msg::kInfo, NuMMParameters::LowerDm2Limit(), NuMMParameters::LowerEpsilonLimit(), and MSG.

Referenced by FCMinimiseNSI().

02960 {
02961   // Calculate the location of the fine grid search
02962   Double_t gridL, gridR, gridT, gridB;
02963   gridL = gridR = gridT = gridB = -1.0;
02964 
02965 //   // Go +- 0.2 in sin space
02966 //   gridL = mmStart.Sn2Bar() - 0.2;
02967 //   gridR = gridL + 0.4;
02968 //   // Now limit it to the physical realm
02969 //   if (gridL <= mmStart.LowerSn2BarLimit()) gridL = mmStart.LowerSn2BarLimit();
02970 //   if (gridR > 1.0) gridR = 1.0;
02971 
02972   // And go +-2e-3 in mass space
02973   gridB = mmStart.Dm2() - 2e-3;
02974   gridT = gridB + 4e-3;
02975   // We only really need to worry about lower limit here
02976   if (gridB <= mmStart.LowerDm2Limit()) gridB = mmStart.LowerDm2Limit();
02977 
02978   // Go +- 0.2 in epsilon space
02979   gridL = mmStart.Epsilon() - 0.2;
02980   gridR = gridL + 0.4;
02981   // Now limit it to the physical realm
02982   if (gridL <= mmStart.LowerEpsilonLimit()) gridL = mmStart.LowerEpsilonLimit();
02983 
02984   MSG("NuSystFitter", Msg::kInfo)
02985     << "Doing fine grid search on eps = [ "
02986     << gridL << ", " << gridR << " ]\n"
02987     << "                              dm  = [ "
02988     << gridB << ", " << gridT << " ]\n";
02989 
02990   // Ensure that we haven't done anything stupid like set low above high
02991   if (gridL >= gridR || gridT <= gridB) {
02992 //      cout << "Error: Negative range in fine grid search" << endl;
02993       throw runtime_error("Negative range in fine grid search");
02994   }
02995 
02996   // Make the surface
02997 
02998   // The width and height of the grid
02999   const Int_t gridWidth = 19;
03000   const Int_t gridHeight = 19;
03001 
03002   // Generate the surface
03003   TH2D likesurf("likesurf", "likesurf", gridWidth, gridL, gridR,
03004     gridHeight, gridB, gridT);
03005 
03006   // And now do a grid search on this
03007   NuMMParameters mmGrid = GridSearchNSI(likesurf, mmStart);
03008 
03009   // Output the minimum
03010   MSG("NusystFitter", Msg::kInfo) << "Fine grid search found epsilon: "
03011     << mmGrid.Epsilon() << ", dm2: " << mmGrid.Dm2() << '\n';
03012 
03013   return mmGrid;
03014 }

NuMMParameters NuSystFitter::FineGridSearchNu ( NuMMParameters  mmStart  )  [protected]

Definition at line 2799 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), MuELoss::e, GridSearchNu(), Msg::kInfo, NuMMParameters::LowerDm2BarLimit(), NuMMParameters::LowerSn2BarLimit(), MSG, NuMMParameters::Sn2(), and NuMMParameters::Sn2Bar().

Referenced by FCMinimiseNu().

02800 {
02801   // Calculate the location of the fine grid search
02802   Double_t gridL, gridR, gridT, gridB;
02803   gridL = gridR = gridT = gridB = -1.0;
02804   
02805   // Go +- 0.2 in sin space
02806   gridL = mmStart.Sn2Bar() - 0.2;
02807   gridR = gridL + 0.4;
02808   // Now limit it to the physical realm
02809   if (gridL <= mmStart.LowerSn2BarLimit()) gridL = mmStart.LowerSn2BarLimit();
02810   if (gridR > 1.0) gridR = 1.0;
02811   
02812   // And go +-2e-3 in mass space
02813   gridB = mmStart.Dm2Bar() - 2e-3;
02814   gridT = gridB + 4e-3;
02815   // We only really need to worry about lower limit here
02816   if (gridB <= mmStart.LowerDm2BarLimit()) gridB = mmStart.LowerSn2BarLimit();
02817   
02818   
02819   MSG("NuSystFitter", Msg::kInfo) 
02820   << "Doing fine grid search on sin = [ "
02821   << gridL << ", " << gridR << " ]\n"
02822   << "                              dm  = [ "
02823   << gridB << ", " << gridT << " ]\n";
02824   
02825   // Ensure that we haven't done anything stupid like set low above high
02826   if (gridL >= gridR || gridT <= gridB) {
02827     //      cout << "Error: Negative range in fine grid search" << endl;
02828     throw runtime_error("Negative range in fine grid search");
02829   }
02830   
02831   // Make the surface
02832   
02833   // The width and height of the grid
02834   const Int_t gridWidth = 19;
02835   const Int_t gridHeight = 19;
02836   
02837   // Generate the surface
02838   TH2D likesurf("likesurf", "likesurf", gridWidth, gridL, gridR,
02839                 gridHeight, gridB, gridT);
02840   
02841   // And now do a grid search on this
02842   NuMMParameters mmGrid = GridSearchNu(likesurf, mmStart);
02843   
02844   // Output the minimum
02845   MSG("NusystFitter", Msg::kInfo) << "Fine grid search found sn2: "
02846   << mmGrid.Sn2() << ", dm2: " << mmGrid.Dm2() << '\n';
02847   
02848   return mmGrid;
02849 }

NuMMParameters NuSystFitter::GridSearch ( TH2D &  likesurf,
NuMMParameters  mmGrid,
NuMMParameters Ref = 0 
) [protected]

Definition at line 2692 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), FillLikelihoodSurfaceBar(), Msg::kDebug, MAXMSG, NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

Referenced by CoarseGridSearch(), and FineGridSearch().

02693 {
02694     Double_t likelihood = FillLikelihoodSurfaceBar(likesurf, mmGrid);
02695 
02696     // Find the point on this surface that is minimum
02697     Int_t minbinX = 0, minbinY = 0, minbinZ = 0;
02698     likesurf.GetMinimumBin(minbinX, minbinY, minbinZ);
02699     Double_t minsin = likesurf.GetXaxis()->GetBinCenter(minbinX);
02700     Double_t minmass = likesurf.GetYaxis()->GetBinCenter(minbinY);
02701 
02702     // Check the zero point, to see if it is lower than anything else on the grid
02703     mmGrid.Sn2Bar(0);
02704     mmGrid.Dm2Bar(0);
02705     Double_t zeropoint = (*this)(mmGrid.VectorParameters());
02706     if (zeropoint < likelihood) {
02707         // Then the no-oscillation case is the minimum
02708         minsin = 0.0;
02709         minmass = 0.0;
02710     }
02711 
02712     // If we have been given a reference point, check to see if this is lower
02713     if (Ref) {
02714       Double_t refpoint = (*this)(Ref->VectorParameters());
02715       if (refpoint < likelihood) {
02716         MAXMSG("NuSystFitter", 5, Msg::kDebug) << "Reference point is lower than grid search" << endl;
02717         minsin = Ref->Sn2Bar();
02718         minmass = Ref->Dm2Bar();
02719       }
02720     }
02721 
02722     // Set mmPars to the newly found point, then return it
02723     mmGrid.Sn2Bar(minsin);
02724     mmGrid.Dm2Bar(minmass);
02725 
02726     return mmGrid;
02727 }

NuMMParameters NuSystFitter::GridSearchNSI ( TH2D &  likesurf,
NuMMParameters  mmGrid,
NuMMParameters Ref = 0 
) [protected]

Definition at line 3017 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Epsilon(), FillLikelihoodSurfaceNSI(), Msg::kDebug, MAXMSG, and NuMMParameters::VectorParameters().

Referenced by CoarseGridSearchNSI(), and FineGridSearchNSI().

03018 {
03019   Double_t likelihood = FillLikelihoodSurfaceNSI(likesurf, mmGrid, 1);
03020 
03021     // Find the point on this surface that is minimum
03022     Int_t minbinX = 0, minbinY = 0, minbinZ = 0;
03023     likesurf.GetMinimumBin(minbinX, minbinY, minbinZ);
03024     Double_t mineps = likesurf.GetXaxis()->GetBinCenter(minbinX);
03025     Double_t minmass = likesurf.GetYaxis()->GetBinCenter(minbinY);
03026 
03027     // Check the zero point, to see if it is lower than anything else on the grid
03028     mmGrid.Epsilon(0);
03029     mmGrid.Dm2(0);
03030     Double_t zeropoint = (*this)(mmGrid.VectorParameters());
03031     if (zeropoint < likelihood) {
03032         // Then the no-oscillation case is the minimum
03033         mineps = 0.0;
03034         minmass = 0.0;
03035     }
03036 
03037     // If we have been given a reference point, check to see if this is lower
03038     if (Ref) {
03039       Double_t refpoint = (*this)(Ref->VectorParameters());
03040       if (refpoint < likelihood) {
03041         MAXMSG("NuSystFitter", 5, Msg::kDebug) << "Reference point is lower than grid search" << endl;
03042         mineps = Ref->Epsilon();
03043         minmass = Ref->Dm2();
03044       }
03045     }
03046 
03047     // Set mmPars to the newly found point, then return it
03048     mmGrid.Epsilon(mineps);
03049     mmGrid.Dm2(minmass);
03050 
03051     return mmGrid;
03052 }

NuMMParameters NuSystFitter::GridSearchNu ( TH2D &  likesurf,
NuMMParameters  mmGrid,
NuMMParameters Ref = 0 
) [protected]

Definition at line 2852 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), FillLikelihoodSurface(), Msg::kDebug, MAXMSG, NuMMParameters::Sn2(), and NuMMParameters::VectorParameters().

Referenced by CoarseGridSearchNu(), and FineGridSearchNu().

02853 {
02854   Double_t likelihood = FillLikelihoodSurface(likesurf, mmGrid);
02855   
02856   // Find the point on this surface that is minimum
02857   Int_t minbinX = 0, minbinY = 0, minbinZ = 0;
02858   likesurf.GetMinimumBin(minbinX, minbinY, minbinZ);
02859   Double_t minsin = likesurf.GetXaxis()->GetBinCenter(minbinX);
02860   Double_t minmass = likesurf.GetYaxis()->GetBinCenter(minbinY);
02861   
02862   // Check the zero point, to see if it is lower than anything else on the grid
02863   mmGrid.Sn2(0);
02864   mmGrid.Dm2(0);
02865   Double_t zeropoint = (*this)(mmGrid.VectorParameters());
02866   if (zeropoint < likelihood) {
02867     // Then the no-oscillation case is the minimum
02868     minsin = 0.0;
02869     minmass = 0.0;
02870   }
02871   
02872   // If we have been given a reference point, check to see if this is lower
02873   if (Ref) {
02874     Double_t refpoint = (*this)(Ref->VectorParameters());
02875     if (refpoint < likelihood) {
02876       MAXMSG("NuSystFitter", 5, Msg::kDebug) << "Reference point is lower than grid search" << endl;
02877       minsin = Ref->Sn2();
02878       minmass = Ref->Dm2();
02879     }
02880   }
02881   
02882   // Set mmPars to the newly found point, then return it
02883   mmGrid.Sn2(minsin);
02884   mmGrid.Dm2(minmass);
02885   
02886   return mmGrid;
02887 }

double NuSystFitter::Likelihood ( const NuMMParameters mmPars  )  const [virtual]

Definition at line 657 of file NuSystFitter.cxx.

References fQuietMode, fRuns, and NuMMParameters::PenaltyTerm().

Referenced by NuABFitter::DeltaChi(), NuABFitter::DoOutput(), NuABFitter::FillLikelihoodSurfaceAB(), FillLikelihoodSurfaceNSIDm2Eps(), FillLikelihoodSurfaceNSIDm2Sn2(), FillLikelihoodSurfaceNSISn2Eps(), FillLikelihoodSurfaceTransCMuMuBoth(), FillLikelihoodSurfaceTransCTauTauBoth(), FillLikelihoodSurfaceTransDm2BarGMu(), FillLikelihoodSurfaceTransDm2BarSn2Bar(), FillLikelihoodSurfaceTransDm2Both(), FillLikelihoodSurfaceTransDm2GMu(), FillLikelihoodSurfaceTransDm2GMuBoth(), FillLikelihoodSurfaceTransDm2NQ(), FillLikelihoodSurfaceTransDm2PQ(), FillLikelihoodSurfaceTransDm2Sn2(), FillLikelihoodSurfaceTransDm2Sn2Both(), FillLikelihoodSurfaceTransGMuBoth(), FillLikelihoodSurfaceTransGMuNQ(), FillLikelihoodSurfaceTransGMuPQ(), FillLikelihoodSurfaceTransGTauBoth(), FillLikelihoodSurfaceTransSn2BarGMu(), FillLikelihoodSurfaceTransSn2Both(), FillLikelihoodSurfaceTransSn2GMu(), FillLikelihoodSurfaceTransSn2GMuBoth(), FillLikelihoodSurfaceTransSn2NQ(), FillLikelihoodSurfaceTransSn2PQ(), NuFCFitterNuMuBar::Fit(), operator()(), NuABFitter::RecoverNegativeDeltaChi(), TransGMuFixDm2Sn2(), and NuFCDelegateArchiver::UseFitter().

00658 {
00659     Double_t likelihood = 0.0;
00660     for (std::vector<NuMMRun*>::iterator runIt = fRuns->begin();
00661          fRuns->end() != runIt;
00662          ++runIt){
00663       likelihood += (*runIt)->ComparePredWithData(mmPars);
00664     }
00665     if (!fQuietMode){
00666       cout << "Stats likelihood: " << likelihood;
00667     }
00668     likelihood += mmPars.PenaltyTerm();
00669     if (!fQuietMode){
00670       cout << " + penalty = " << likelihood << endl;
00671     }
00672 
00673     return likelihood;
00674 }

TH2D * NuSystFitter::LikelihoodSurfaceBar ( const NuMMParameters mmPars,
Int_t  sinPoints,
Int_t  massPoints,
Double_t  Sn2BarMin = 0.0,
Double_t  Sn2BarMax = 1.0,
Double_t  Dm2BarMin = 0.0,
Double_t  Dm2BarMax = 12e-3,
Bool_t  True_Minimum = true 
) [virtual]

Function to generate a likelihood surface.

Definition at line 1116 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), FillLikelihoodSurfaceBar(), fQuietMode, Minimise(), minimum, NuMMParameters::MinuitPass(), NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

01124 {
01125     // Calculate the bin spacing. This is -1 because we want to 'add'
01126     // one bin outside of the specified range, so that the bin contents
01127     // are representative of the center of the bin.
01128     Double_t sininterval =  (sn2BarMax - sn2BarMin) / (sinPoints -1);
01129     Double_t massinterval = (dm2BarMax - dm2BarMin) / (massPoints - 1);
01130 
01131     // Offsets. This is how much to shift the bins min and max by
01132     Double_t sinoffset = sininterval / 2.0;
01133     Double_t massoffset = massinterval / 2.0;
01134     // Create the destination histogram
01135     TH2D *surface = new TH2D(   "chisurf", "#chi^{2} Likelihood Surface",
01136                                 sinPoints, sn2BarMin-sinoffset, sn2BarMax+sinoffset,
01137                                 massPoints, dm2BarMin-massoffset, dm2BarMax+massoffset );
01138     // Set it up
01139     surface->GetXaxis()->SetTitle("sin^{2}2#bar#theta");
01140     surface->GetYaxis()->SetTitle("#Delta#barm^{2}");
01141 
01142     // A variable to store the minimum likelihood encountered
01143     // Double_t minimum = -1;
01144 
01145     // Make a copy of the parameters
01146     NuMMParameters mmPars(mmParsIn);
01147 
01148     // Variables to store the current value. These are declared here
01149     // so that we can verify we have counted them correctly after the
01150     // looping ends.
01151     // Double_t sinVal = 0., dmVal = 0.;
01152 
01153     if (!fQuietMode) cout << "Generating Likelihood Surface..." << endl;
01154 
01155     Double_t minimum = this->FillLikelihoodSurfaceBar(*surface, mmPars);
01156 
01157     // // Now, start looping over the surface
01158     // for (Int_t sinBin = 1; sinBin <= sinPoints; sinBin++)
01159     // {
01160     //     // Calculate and set this sin (bin 0 starts at -sinoffset but
01161     //     // we want it to contain the value at 0)
01162     //     sinVal = sn2BarMin + sininterval*(sinBin-1);
01163     //     mmPars.Sn2Bar(sinVal);
01164     //
01165     //     for (Int_t massBin = 1; massBin <= massPoints; massBin++)
01166     //     {
01167     //         // Calculate and set the dmass value
01168     //         dmVal = dm2BarMin + massinterval*(massBin-1);
01169     //         mmPars.Dm2Bar(dmVal);
01170     //
01171     //         // Now, calculate the likelihood using this fitter
01172     //         Double_t likelihood = (*this)(mmPars.VectorParameters());
01173     //
01174     //         // Minimim calculations
01175     //         if (minimum < 0 || likelihood < minimum)
01176     //             minimum = likelihood;
01177     //
01178     //         // Set the bin content!
01179     //         surface->SetBinContent(sinBin, massBin, likelihood);
01180     //     }
01181     // }
01182     // Print a message to say we have finished
01183     if (!fQuietMode) {
01184         cout << "Surface Generation Done. Minimum Likelihood = "
01185              << minimum << endl;
01186     }
01187 
01188 
01189     // Get the minimum bin on this histogram
01190     Int_t minbinX, minbinY, minbinZ;
01191     surface->GetMinimumBin(minbinX, minbinY, minbinZ);
01192     // Set the fit parameters to this point as a start point
01193     mmPars.Sn2Bar(surface->GetXaxis()->GetBinCenter(minbinX));
01194     mmPars.Dm2Bar(surface->GetYaxis()->GetBinCenter(minbinY));
01195 
01196   if (True_Minimum) {
01197     // Now minimise at this point, and grab the likelihood
01198     if (!fQuietMode) cout << "Searching for true minimum" << endl;
01199 
01200     NuMMParameters mmFit = this->Minimise(mmPars);
01201     Double_t minlike = (*this)(mmFit.VectorParameters());
01202 
01203     // Check if the fit failed, if it did, then don't use it
01204     if (mmFit.MinuitPass())
01205     {
01206         // If this is smaller than this minimum we have found (it should be)
01207         if (minlike <= minimum)
01208           minimum = minlike;
01209         else
01210           cout << "WARNING: Fit Minimisation found minimum higher than grid search" << endl;
01211     }
01212   }
01213     // Now, we have generated likelihoods for all grid points.
01214     // Offset each one by the minimum, so that it is a delta chi
01215     // squared surface.
01216     for (Int_t sinBin = 1; sinBin <= sinPoints; sinBin++) {
01217         for (Int_t massBin = 1; massBin <= massPoints; massBin++) {
01218             Double_t content = surface->GetBinContent(sinBin, massBin);
01219             surface->SetBinContent(sinBin, massBin, content - minimum);
01220         }
01221     }
01222 
01223     // Now we have generated and offset the histogram. Return it!
01224     return surface;
01225 }

const NuMMParameters NuSystFitter::Minimise ( const NuMMParameters mmPars,
bool &  success 
) [virtual]

Definition at line 104 of file NuSystFitter.cxx.

References NuMMParameters::AreAllParametersFixed(), NuMMParameters::CopyConstraintsLimitsAndFixings(), fQuietMode, NuMMParameters::kAlpha, NuMMParameters::kDecayPipe, NuMMParameters::kDelta1, NuMMParameters::kDelta2, NuMMParameters::kDelta3, NuMMParameters::kDm2, NuMMParameters::kDm221, NuMMParameters::kDm243, NuMMParameters::kDm2Bar, NuMMParameters::kEpsilon, NuMMParameters::kgmu, NuMMParameters::kMu2, NuMMParameters::kNcBack, NuMMParameters::kNorm, NuMMParameters::kNormTrans, NuMMParameters::kNuBarXSecDIS, NuMMParameters::kNumParameters, NuMMParameters::kShwScale, NuMMParameters::kShwScaleTrans, NuMMParameters::kSn2, NuMMParameters::kSn2Bar, NuMMParameters::kTheta12, NuMMParameters::kTheta13, NuMMParameters::kTheta14, NuMMParameters::kTheta24, NuMMParameters::kTheta34, NuMMParameters::kTrkScaleTrans, NuMMParameters::MinuitParameters(), NuMMParameters::MinuitPass(), and NuMMParameters::PrintStatus().

00106 {
00107   success = true;
00108 
00109   if(mmPars.AreAllParametersFixed()) return mmPars;
00110 
00111   //Set up Minuit parameters
00112   ROOT::Minuit2::MnUserParameters mnPars = mmPars.MinuitParameters();
00113   std::cout << "eps: " << mnPars.Precision().Eps() << std::endl;; //EPS2()
00114 //  mnPars.SetPrecision(1e-11);
00115 
00116   vector<double> param_vec = mnPars.Params();
00117 //    cout << "Fitting with " << param_vec.size() << " parameters." << endl;
00118 //   cout << "Fitting with " << mnPars.VariableParameters() << " variable parameters." << endl;
00119   if (!fQuietMode) mmPars.PrintStatus();
00120 
00121   //Run Minuit fit
00122   ROOT::Minuit2::MnMigrad mnMigrad(*this,mnPars);
00123   ROOT::Minuit2::FunctionMinimum funcMin = mnMigrad();
00124 //  ROOT::Minuit2::FunctionMinimum funcMin = mnMigrad(0,1);
00125 //  ROOT::Minuit2::MnSimplex mnSimplex(*this,mnPars);
00126 //  ROOT::Minuit2::FunctionMinimum funcMin = mnSimplex();
00127 //  ROOT::Minuit2::MnMinimize mnMinimize(*this,mnPars);
00128 //  ROOT::Minuit2::FunctionMinimum funcMin = mnMinimize();
00129 
00130   if (!funcMin.IsValid()) {
00131     success = false;
00132     cout << "Function minimum is not valid!" << endl;
00133   } else {
00134     success = true;
00135     if (!fQuietMode){
00136       const ROOT::Minuit2::MnUserParameters upars = funcMin.UserParameters();
00137           const ROOT::Minuit2::MnUserCovariance covariance = funcMin.UserCovariance();
00138 
00139           cout << funcMin << endl;
00140           cout << covariance << endl;
00141 
00142       cout << "Function minimum is: " << endl
00143            << "Sn2: "       << upars.Parameter(NuMMParameters::kSn2).Value()
00144            << ", dm2: "     << upars.Parameter(NuMMParameters::kDm2).Value()
00145            << "; norm: "    << upars.Parameter(NuMMParameters::kNorm).Value()
00146            << "; shwEn: "   << upars.Parameter(NuMMParameters::kShwScale).Value()
00147            << "; NCBack: "  << upars.Parameter(NuMMParameters::kNcBack).Value()
00148            << "; Sn2bar: "  << upars.Parameter(NuMMParameters::kSn2Bar).Value()
00149            << "; dm2bar: "  << upars.Parameter(NuMMParameters::kDm2Bar).Value()
00150            << "; epsilon: " << upars.Parameter(NuMMParameters::kEpsilon).Value()
00151            << "; Theta12: " << upars.Parameter(NuMMParameters::kTheta12).Value()
00152            << "; Theta13: " << upars.Parameter(NuMMParameters::kTheta13).Value()
00153            << "; Theta14: " << upars.Parameter(NuMMParameters::kTheta14).Value()
00154            << "; Delta1: "  << upars.Parameter(NuMMParameters::kDelta1).Value()
00155            << "; Delta2: "  << upars.Parameter(NuMMParameters::kDelta2).Value()
00156            << "; Delta3: "  << upars.Parameter(NuMMParameters::kDelta3).Value()
00157            << "; Theta24: " << upars.Parameter(NuMMParameters::kTheta24).Value()
00158            << "; Theta34: " << upars.Parameter(NuMMParameters::kTheta34).Value()
00159            << "; Dm221: "   << upars.Parameter(NuMMParameters::kDm221).Value()
00160            << "; Dm243: "   << upars.Parameter(NuMMParameters::kDm243).Value()
00161            << "; mu2: "     << upars.Parameter(NuMMParameters::kMu2).Value()
00162            << "; alpha: "   << upars.Parameter(NuMMParameters::kAlpha).Value()
00163            << "; CPT: "     << upars.Parameter(NuMMParameters::kNumParameters).Value()
00164            << "; normTrans: "    << upars.Parameter(NuMMParameters::kNormTrans).Value()
00165            << "; shwEnTrans: "   << upars.Parameter(NuMMParameters::kShwScaleTrans).Value()
00166            << "; trkEnTrans: "   << upars.Parameter(NuMMParameters::kTrkScaleTrans).Value()
00167            << "; decayPipe: "    << upars.Parameter(NuMMParameters::kDecayPipe).Value()
00168            << "; nubarXSecDIS: " << upars.Parameter(NuMMParameters::kNuBarXSecDIS).Value()
00169            << "; gmu: "          << upars.Parameter(NuMMParameters::kgmu).Value()
00170       << endl;
00171     }
00172   }
00173 
00174   //cout << "Diagonal covar elements:" << endl;
00175   //ROOT::Minuit2::MnUserCovariance covar = funcMin.UserCovariance();
00176   //for (unsigned int i = 0; i < covar.Nrow(); i++) {
00177   //  cout << "element " << i << " = " << covar(i,i) << endl;
00178   //}
00179 
00180   NuMMParameters bestFitPars(funcMin);
00181   // Make sure the parameters we return are configured the same way as what was passed
00182   bestFitPars.CopyConstraintsLimitsAndFixings(mmPars);
00183   bestFitPars.MinuitPass(success);
00184   //bestFitPars.PrintStatus();
00185   return bestFitPars;
00186 }

const NuMMParameters NuSystFitter::Minimise ( const NuMMParameters mmPars  )  [virtual]

Definition at line 97 of file NuSystFitter.cxx.

Referenced by FeldmanSterile::CalculateChi2Best2D(), FeldmanSterile::CalculateChi2Profile(), Chi2ValleyBar(), Chi2ValleyBarDm2Bar(), Chi2ValleyDm2(), Chi2ValleySn2(), FarOverNearFitPQ::DoFit(), FarOverNearFit::DoFit(), FCMinimise(), FCMinimiseNSI(), FCMinimiseNu(), FillCPT(), FillLikelihoodSurfaceGeneric(), FillLikelihoodSurfaceNSIDm2Eps(), FillLikelihoodSurfaceNSIDm2Sn2(), FillLikelihoodSurfaceNSISn2Eps(), FillLikelihoodSurfaceTransCMuMuBoth(), FillLikelihoodSurfaceTransCTauTauBoth(), FillLikelihoodSurfaceTransDm2BarGMu(), FillLikelihoodSurfaceTransDm2BarSn2Bar(), FillLikelihoodSurfaceTransDm2Both(), FillLikelihoodSurfaceTransDm2GMu(), FillLikelihoodSurfaceTransDm2GMuBoth(), FillLikelihoodSurfaceTransDm2NQ(), FillLikelihoodSurfaceTransDm2PQ(), FillLikelihoodSurfaceTransDm2Sn2(), FillLikelihoodSurfaceTransDm2Sn2Both(), FillLikelihoodSurfaceTransGMuBoth(), FillLikelihoodSurfaceTransGMuNQ(), FillLikelihoodSurfaceTransGMuPQ(), FillLikelihoodSurfaceTransGTauBoth(), FillLikelihoodSurfaceTransSn2BarGMu(), FillLikelihoodSurfaceTransSn2Both(), FillLikelihoodSurfaceTransSn2GMu(), FillLikelihoodSurfaceTransSn2GMuBoth(), FillLikelihoodSurfaceTransSn2NQ(), FillLikelihoodSurfaceTransSn2PQ(), NuFCFitterNuMuBar::Fit(), NuEZFitterNSI::Fit(), NuEZFitter::Fit(), FarOverNearFitPQ::FitTheta34Only(), FarOverNearFit::FitTheta34Only(), FeldmanSterile::FitThreeFlavourOnly(), LikelihoodSurfaceBar(), NuFCFitterNSINubar::MinuitFit(), NuFCFitterNSINu::MinuitFit(), NuFCFitterNSI::MinuitFit(), NuFCFitter::MinuitFit(), FarOverNearFitPQ::PerformInitialThreeFlavourFit(), FarOverNearFit::PerformInitialThreeFlavourFit(), NuFCFitterNSINubar::RecoverNegativeDeltaChi(), NuFCFitterNSINu::RecoverNegativeDeltaChi(), NuFCFitterNSI::RecoverNegativeDeltaChi(), NuFCFitter::RecoverNegativeDeltaChi(), and TransGMuFixDm2Sn2().

00098 {
00099   bool success;
00100   return Minimise(mmPars, success);
00101 }

const std::vector< std::pair< Double_t, Double_t > > NuSystFitter::NeutrinoContourFinder ( const NuMMParameters mmPars  )  [virtual]

Definition at line 323 of file NuSystFitter.cxx.

References NuMMParameters::FixSn2(), NuMMParameters::Sn2(), and NuMMParameters::VectorParameters().

00324 {
00325   //Find the target chi2 for the contour
00326   NuMMParameters bestFitPars = this->Minimise(inputPars);
00327   Double_t bestFitChi2 = (*this)(bestFitPars.VectorParameters());
00328   Double_t targetChi2 = bestFitChi2 + this->Up();
00329   cout << "Target chi2 = " << targetChi2 << endl;
00330 
00331   //Fix sn2 ready for the dm2 scans:
00332   NuMMParameters mmPars = inputPars;
00333   mmPars.FixSn2();
00334   //Make somewhere to put the contour points
00335   std::vector<pair<Double_t,Double_t> > pointsUp;
00336   std::vector<pair<Double_t,Double_t> > pointsDown;
00337 
00338   //Run the contour in quiet mode
00339   this->QuietModeOn();
00340 
00341   //Get the points for the contour
00342   for (vector<Double_t>::iterator snIt = fsn2PointsForNeutrinoContour.begin();
00343        snIt != fsn2PointsForNeutrinoContour.end();
00344        ++snIt){
00345     Double_t sn2 = *snIt;
00346     mmPars.Sn2(sn2);
00347     const std::pair<Double_t, Double_t> values =
00348       this->DmScanForContour(mmPars,targetChi2);
00349     cout << "Contour points: sn2 = " << sn2
00350          << "; dm2 = " << values.first << " and " << values.second
00351          << endl;
00352     if ((values.first < 0.0) || (values.second < 0.0)){
00353       cout << "End of contour at sn2 = " << mmPars.Sn2() << endl;
00354       break;
00355     }
00356     else{
00357       pair<Double_t, Double_t> pointDown(sn2,values.first);
00358       pair<Double_t, Double_t> pointUp(sn2,values.second);
00359       pointsDown.push_back(pointDown);
00360       pointsUp.push_back(pointUp);
00361     }
00362   }
00363 
00364   //Make a single vector of points
00365   std::vector<pair<Double_t,Double_t> > contour;
00366   for (vector<pair<Double_t, Double_t> >::const_iterator it = pointsUp.begin();
00367        it != pointsUp.end();
00368        ++it){
00369     pair<Double_t,Double_t> currentPoint((*it).first,(*it).second);
00370     contour.push_back(currentPoint);
00371   }
00372   for (vector<pair<Double_t, Double_t> >::reverse_iterator
00373          itr = pointsDown.rbegin();
00374        itr != pointsDown.rend();
00375        ++itr){
00376     pair<Double_t,Double_t> currentPoint((*itr).first,(*itr).second);
00377     contour.push_back(currentPoint);
00378   }
00379   return contour;
00380 }

void NuSystFitter::NoOscPrediction (  )  [virtual]

Definition at line 91 of file NuSystFitter.cxx.

00092 {
00093 
00094 }

virtual void NuSystFitter::NumberOfContourPoints ( const Int_t  numContourPoints  )  [inline, virtual]

Definition at line 47 of file NuSystFitter.h.

References fNumContourPoints.

00048     {fNumContourPoints = numContourPoints;}

const TGraph * NuSystFitter::OneSidedNeutrinoContourFinder ( const NuMMParameters inputPars  )  [virtual]

Definition at line 252 of file NuSystFitter.cxx.

00253 {
00254   std::vector<std::pair<Double_t,Double_t> > contPoints =
00255     this->NeutrinoContourFinder(inputPars);
00256 
00257   //Put the contour points in a TGraph
00258   TGraph* contour = new TGraph(contPoints.size());
00259   Int_t point=0;
00260   for (vector<pair<Double_t, Double_t> >::const_iterator it =
00261          contPoints.begin();
00262        it != contPoints.end();
00263        ++it){
00264     contour->SetPoint(point,(*it).first,(*it).second);
00265     cout << "Set contour point " << point
00266          << ", " << (*it).first
00267          << ", " << (*it).second
00268          << endl;
00269     ++point;
00270   }
00271   return contour;
00272 }

double NuSystFitter::operator() ( const std::vector< double > &  pars  )  const [virtual]

Definition at line 648 of file NuSystFitter.cxx.

References Likelihood().

00649 {
00650   NuMMParameters mmPars(pars);
00651   //cout << "Called operator()" << endl;
00652   //mmPars.PrintStatus();
00653   return Likelihood(mmPars);
00654 }

const TGraph * NuSystFitter::PlotContours ( const NuMMParameters mmPars  )  [virtual]

Definition at line 189 of file NuSystFitter.cxx.

References NuMMParameters::ContourPars(), fNumContourPoints, NuMMParameters::kAlpha, NuMMParameters::kDm2, NuMMParameters::kDm2Bar, NuMMParameters::kMu2, NuMMParameters::kNcBack, NuMMParameters::kNorm, NuMMParameters::kNumParameters, NuMMParameters::kShwScale, NuMMParameters::kSn2, NuMMParameters::kSn2Bar, NuMMParameters::MinuitParameters(), and Munits::second.

00190 {
00191 
00192 
00193   //Set up Minuit parameters
00194   ROOT::Minuit2::MnUserParameters mnPars = mmPars.MinuitParameters();
00195 
00196   //Run Minuit fit
00197   ROOT::Minuit2::MnMigrad mnMigrad(*this,mnPars);
00198   ROOT::Minuit2::FunctionMinimum funcMin = mnMigrad();
00199 
00200   if (!funcMin.IsValid()){
00201     cout << "Function minimum is not valid!" << endl;
00202   }
00203   else{
00204     cout << "Contour Minimization: " << endl;
00205 
00206     const ROOT::Minuit2::MnUserParameters upars = funcMin.UserParameters();
00207 
00208     cout << "Function minimum is: " << endl
00209          << "Sn2: " << upars.Parameter(NuMMParameters::kSn2).Value()
00210          << ", dm2: " << upars.Parameter(NuMMParameters::kDm2).Value()
00211          << "; norm: " << upars.Parameter(NuMMParameters::kNorm).Value()
00212          << "; shwEn: " << upars.Parameter(NuMMParameters::kShwScale).Value()
00213          << "; NCBack: " << upars.Parameter(NuMMParameters::kNcBack).Value()
00214          << "; Sn2bar: " << upars.Parameter(NuMMParameters::kSn2Bar).Value()
00215          << "; dm2bar: " << upars.Parameter(NuMMParameters::kDm2Bar).Value()
00216          << "; mu2: " << upars.Parameter(NuMMParameters::kMu2).Value()
00217          << "; alpha: " << upars.Parameter(NuMMParameters::kAlpha).Value()
00218          << "; CPT: " << upars.Parameter(NuMMParameters::kNumParameters).Value()
00219          << endl;
00220   }
00221 
00222   ROOT::Minuit2::MnContours contours(*this,funcMin);
00223   ROOT::Minuit2::ContoursError contErr = contours.Contour(mmPars.ContourPars().first,
00224                                                           mmPars.ContourPars().second,
00225                                                           fNumContourPoints);
00226 
00227   std::vector<std::pair<Double_t,Double_t> > vContPoints = contErr();
00228   //std::vector<std::pair<Double_t,Double_t> > vContPoints(10);// = contErr();
00229 
00230   Double_t* x = new Double_t[vContPoints.size()];
00231   for (UInt_t pointCount = 0; pointCount<vContPoints.size(); ++pointCount){
00232     x[pointCount] = vContPoints[pointCount].first;
00233   }
00234 
00235   Double_t* y = new Double_t[vContPoints.size()];
00236   for (UInt_t pointCount = 0; pointCount<vContPoints.size(); ++pointCount){
00237     y[pointCount] = vContPoints[pointCount].second;
00238   }
00239   TGraph* gCont = new TGraph(vContPoints.size(),y,x);
00240 //  cout << "drawing contour" << endl;
00241 //  gCont->SetName("gCont");
00242 //  gDirectory->Print();
00243 //  TFile* file=TFile::Open("contourplot.root","RECREATE");
00244 //  gDirectory->Print();
00245 //  gCont->Write("gCont");
00246 //  file->Close();
00247   return gCont;
00248 }

void NuSystFitter::push_back ( NuMMRun run  )  [virtual]

Definition at line 66 of file NuSystFitter.cxx.

References fRuns, and run().

Referenced by NuFCFittingInterface::AddRun(), NuEZRunsFitter::BuildFitter(), FeldmanSterile::Calculate2DdeltaChi2(), NuFCFitter::NuFCFitter(), NuFCFitterNSI::NuFCFitterNSI(), NuFCFitterNSINu::NuFCFitterNSINu(), NuFCFitterNSINubar::NuFCFitterNSINubar(), NuEZFitterNSI::Rebuild(), and NuEZFitter::Rebuild().

00067 {
00068   fRuns->push_back(run);
00069 }

void NuSystFitter::QuietModeOff (  )  [virtual]

Definition at line 695 of file NuSystFitter.cxx.

References fQuietMode, and fRuns.

00696 {
00697   fQuietMode = 0;
00698   for (std::vector<NuMMRun*>::iterator runIt = fRuns->begin();
00699        fRuns->end() != runIt;
00700        ++runIt){
00701     (*runIt)->QuietModeOff();
00702   }
00703 }

void NuSystFitter::QuietModeOn (  )  [virtual]

Definition at line 684 of file NuSystFitter.cxx.

References fQuietMode, and fRuns.

Referenced by BatchModeOn(), and FeldmanSterile::Calculate2DdeltaChi2().

00685 {
00686   fQuietMode = 1;
00687   for (std::vector<NuMMRun*>::iterator runIt = fRuns->begin();
00688        fRuns->end() != runIt;
00689        ++runIt){
00690     (*runIt)->QuietModeOn();
00691   }
00692 }

virtual void NuSystFitter::SetErrorDef ( const Double_t  up  )  [inline, virtual]

Definition at line 25 of file NuSystFitter.h.

References fUp.

00025 {fUp = up;}

virtual void NuSystFitter::Setup ( const std::vector< std::string >  helperInputs,
const std::vector< NuHistInterpolator * >  fdDataInput,
const std::vector< NuHistInterpolator * >  ndDataInput,
const std::string  xsecFile,
const std::vector< std::string >  ,
const std::string  ,
Int_t   
) [virtual]

const std::pair< Double_t, Double_t > NuSystFitter::SingleParDm2BarErrors ( const NuMMParameters mmPars,
const bool  quiet = false 
) [virtual]

Definition at line 503 of file NuSystFitter.cxx.

00504 {
00505   return SingleParErrors(mmPars,5,quiet);
00506 }

const std::pair< Double_t, Double_t > NuSystFitter::SingleParDm2Errors ( const NuMMParameters mmPars,
const bool  quiet = false 
) [virtual]

Definition at line 496 of file NuSystFitter.cxx.

00497 {
00498   return SingleParErrors(mmPars,0,quiet);
00499 }

const std::pair< Double_t, Double_t > NuSystFitter::SingleParErrors ( const NuMMParameters mmPars,
const int  idx,
const bool  quiet 
) [virtual]

Definition at line 517 of file NuSystFitter.cxx.

References NuMMParameters::AreAllParametersFixed(), MuELoss::e, NuMMParameters::FixParameterByIndex(), NuMMParameters::GetParameterByIndex(), NuMMParameters::MinuitPass(), NuMMParameters::SetParameterByIndex(), and NuMMParameters::VectorParameters().

00519 {
00520   std::pair<Double_t,Double_t> errors;
00521 
00522   NuMMParameters bestFitPoint = this->Minimise(mmPars);
00523   Double_t bestFitChi2 = (*this)(bestFitPoint.VectorParameters());
00524   Double_t targetChi2 = bestFitChi2 + this->Up();
00525   Double_t bestFitPar = bestFitPoint.GetParameterByIndex(idx);
00526 
00527   cout << "Best chi2: " << bestFitChi2
00528        << " at par = " << bestFitPar
00529        << ". Target chi2: " << targetChi2
00530        << endl;
00531 
00532   for (Int_t direction = -1; direction < 2; direction += 2){
00533 
00534     if(!quiet) cout << "Direction = " << direction << endl;
00535 
00536     NuMMParameters scanMMPars;
00537 
00538 //    Double_t incrementArray[5] = {0.1e-3,
00539 //                                0.01e-3,
00540 //                                0.001e-3,
00541 //                                0.0001e-3,
00542 //                                0.00001e-3};
00543     Double_t firstChi2 = bestFitChi2;
00544     Double_t firstPar = bestFitPar;
00545     Double_t currentChi2 = -1.0;
00546     Double_t lastChi2 = -1.0;
00547     Double_t lastPar = -1.0;
00548     Double_t increment = 1e-3*bestFitPar;
00549 
00550     while (currentChi2<targetChi2 && increment < 1e6){
00551       increment *= 2.0;
00552       scanMMPars = mmPars;
00553       scanMMPars.FixParameterByIndex(idx);
00554       scanMMPars.SetParameterByIndex(idx,bestFitPar + direction*increment);
00555       lastPar = scanMMPars.GetParameterByIndex(idx);
00556       if(!quiet) cout << "Trying par = " << scanMMPars.GetParameterByIndex(idx) << "... ";
00557       if (scanMMPars.AreAllParametersFixed()){
00558         currentChi2 = (*this)(scanMMPars.VectorParameters());
00559       }
00560       else{
00561         currentChi2 = (*this)(this->Minimise(scanMMPars).VectorParameters());
00562       }
00563       if(!quiet){
00564         cout << "Starting while with increment = " << increment
00565              << ", currentChi2 = " << currentChi2
00566              << ", lastPar = " << lastPar
00567              << endl;
00568       }
00569     }
00570 
00571     increment /= 2.0;
00572 
00573     while (TMath::Abs(currentChi2-targetChi2)>1e-5 && increment > 1e-5*bestFitPar){
00574 
00575       scanMMPars = mmPars;
00576       scanMMPars.FixParameterByIndex(idx);
00577       scanMMPars.SetParameterByIndex(idx,firstPar);
00578 
00579       if(!quiet){
00580         cout << "Starting while with firstChi2 = " << firstChi2
00581              << ", firstPar = " << firstPar
00582              << ", increment = " << increment
00583              << endl;
00584       }
00585 
00586       lastChi2 = currentChi2;
00587       lastPar = scanMMPars.GetParameterByIndex(idx);
00588       Double_t passminuit = true;
00589       scanMMPars.SetParameterByIndex(idx,scanMMPars.GetParameterByIndex(idx) + direction*increment);
00590       if(!quiet) cout << "Trying par = " << scanMMPars.GetParameterByIndex(idx) << "... ";
00591       if (scanMMPars.AreAllParametersFixed()){
00592         currentChi2 = (*this)(scanMMPars.VectorParameters());
00593       }
00594       else{
00595         NuMMParameters minPars = Minimise(scanMMPars);
00596         currentChi2 = (*this)(minPars.VectorParameters());
00597         passminuit = minPars.MinuitPass();
00598       }
00599       if(!quiet) cout << "chi2 = " << currentChi2 << endl;
00600 
00601       if(passminuit){
00602         if(currentChi2 > targetChi2) firstPar = lastPar;
00603         else firstPar = lastPar + direction*increment;
00604         increment /= 2.0;
00605       }
00606       else{
00607         firstPar = lastPar;
00608         increment *= 1.01;
00609       }
00610     }
00611 
00612     Double_t chi2Difference = lastChi2 - currentChi2;
00613     Double_t parDifference = lastPar - scanMMPars.GetParameterByIndex(idx);
00614     Double_t value = -1.0;
00615     if(!quiet){
00616       cout << "chi2Difference = " << -chi2Difference << endl;
00617       cout << "parDifference = " << parDifference << endl;
00618     }
00619     if ((fabs(chi2Difference)<1e-12) || (fabs(parDifference)<1e-12)){
00620       cout << "Chi2 or Par difference is zero. Supplying currentChi2."
00621            << endl;
00622       value = scanMMPars.GetParameterByIndex(idx);
00623     }
00624     else{
00625       Double_t fractionBetween = (targetChi2 - currentChi2)/chi2Difference;
00626       value = fractionBetween*parDifference + scanMMPars.GetParameterByIndex(idx);
00627     }
00628     if(!quiet){
00629       cout << "Error is between " << lastPar
00630            << " and " << scanMMPars.GetParameterByIndex(idx)
00631            << ", with chi2 between " << lastChi2
00632            << " and " << currentChi2
00633            << endl;
00634       cout << "Value is " << value << endl;
00635     }
00636     if (-1 == direction){errors.first = value;}
00637     if (1 == direction){errors.second = value;}
00638   }
00639 
00640   cout << "Final answer is par = " << bestFitPar
00641        << " + " << errors.second - bestFitPar
00642        << " - " << -errors.first + bestFitPar
00643        << endl;
00644   return errors;
00645 }

const std::pair< Double_t, Double_t > NuSystFitter::SingleParSn2BarErrors ( const NuMMParameters mmPars,
const bool  quiet = false 
) [virtual]

Definition at line 510 of file NuSystFitter.cxx.

00511 {
00512   return SingleParErrors(mmPars,6,quiet);
00513 }

NuMMParameters NuSystFitter::Sn2GridSearch ( Double_t  Dm2,
Double_t  Dm2Bar,
Int_t  resolution = 0 
) [protected]

Does a grid search over sn2/sn2bar, with resolution bins per parameter.

Definition at line 3248 of file NuSystFitter.cxx.

References NuMMParameters::Dm2(), NuMMParameters::Dm2Bar(), NuMMParameters::Sn2(), NuMMParameters::Sn2Bar(), and NuMMParameters::VectorParameters().

Referenced by FillCPT().

03249 {
03250   if (resolution <= 0) resolution = 11;
03251   // Use a 2D histogram to define the search grid
03252   TH2D sinsurf("sinsurf","sinsurf",resolution,0,1,resolution,0,1);
03253   // Create a template parameters object
03254   NuMMParameters pars;
03255   pars.Dm2(Dm2);
03256   pars.Dm2Bar(Dm2Bar);
03257 
03258   Double_t minLike, minSn2, minSn2Bar;
03259   minLike = minSn2 = minSn2Bar = -1;
03260 
03261   // Loop over every bin
03262   for (Int_t xBin = 1; xBin <= sinsurf.GetNbinsX(); ++xBin) {
03263     for (Int_t yBin = 1; yBin <= sinsurf.GetNbinsY(); ++yBin) {
03264       Double_t sn2 = sinsurf.GetYaxis()->GetBinCenter(yBin);
03265       Double_t sn2bar = sinsurf.GetXaxis()->GetBinCenter(xBin);
03266       pars.Sn2(sn2);
03267       pars.Sn2Bar(sn2bar);
03268 
03269       // Evaluate the likelihood
03270       const Double_t likelihood = (*this)(pars.VectorParameters());
03271 
03272       // Check for minimum (or first entry)
03273       if ((xBin == 1 && yBin == 1) || likelihood < minLike) {
03274         minLike = likelihood;
03275         minSn2 = sn2;
03276         minSn2Bar = sn2bar;
03277       }
03278     }
03279   }
03280 
03281   // Check Zero and One
03282   pars.Sn2(0.0);
03283   pars.Sn2Bar(0.0);
03284   Double_t likelihood = (*this)(pars.VectorParameters());
03285   if (likelihood < minLike) {
03286     minLike = likelihood;
03287     minSn2 = 0.0;
03288     minSn2Bar = 0.0;
03289   }
03290 
03291   // Check Zero and One
03292   pars.Sn2(1.0);
03293   pars.Sn2Bar(1.0);
03294   likelihood = (*this)(pars.VectorParameters());
03295   if (likelihood < minLike) {
03296     minLike = likelihood;
03297     minSn2 = 1.0;
03298     minSn2Bar = 1.0;
03299   }
03300 
03301   // We should now have a minimum. Sanity check this
03302   assert(minLike > 0);
03303   assert(minSn2 > 0);
03304   assert(minSn2Bar > 0);
03305 
03306   // Now build our return
03307   pars.Sn2(minSn2);
03308   pars.Sn2Bar(minSn2Bar);
03309 
03310   return pars;
03311 }

virtual void NuSystFitter::Sn2PointsForNeutrinoContour ( const std::vector< Double_t > &  sn2PointsForNeutrinoContour  )  [inline, virtual]

Definition at line 62 of file NuSystFitter.h.

References fsn2PointsForNeutrinoContour.

00063     {fsn2PointsForNeutrinoContour = sn2PointsForNeutrinoContour;}

Double_t NuSystFitter::TransGMuFixDm2Sn2 ( TH1D &  surface,
TH1D &  dm2vals,
TH1D &  sn2vals,
NuMMParameters  mmPars 
) [virtual]

Definition at line 2266 of file NuSystFitter.cxx.

References NuMMParameters::Dm2Bar(), NuMMParameters::FixGMu(), NuMMParameters::GMu(), Likelihood(), Minimise(), minimum, and NuMMParameters::Sn2Bar().

02270 {
02271 
02272   Int_t Ngmu = surface.GetNbinsX();
02273 
02274   Double_t minimum = 1e6;
02275 
02276   for (int i = 1; i <= Ngmu; i++){
02277 
02278     //  if(fQuietMode < 2){
02279     //    NuUtilities::ProgressBar((i-1)*Ndm2+j, Ndm2*Ngmu, 1);
02280     //  }
02281 
02282       Double_t gmu = surface.GetXaxis()->GetBinCenter(i);
02283       mmPars.GMu(gmu);
02284       mmPars.FixGMu();
02285 
02286 //      mmPars.FixDm2Bar();
02287 //      mmPars.FixSn2Bar();
02288 //      mmPars.FixDm2();
02289 //      mmPars.FixSn2();
02290 //      mmPars.ConstrainCPT();
02291 
02292       NuMMParameters bf = this->Minimise(mmPars);
02293       Double_t likelihood = Likelihood(bf);
02294       cout<<"BestFit Likelihood: ="<<likelihood<<endl;
02295       Double_t sn2Min = bf.Sn2Bar();
02296       Double_t dm2Min = bf.Dm2Bar();
02297 
02298       surface.SetBinContent(i, likelihood);
02299       dm2vals.SetBinContent(i, dm2Min);
02300       sn2vals.SetBinContent(i, sn2Min);
02301       if(likelihood < minimum) minimum = likelihood;
02302   }
02303 
02304   return minimum;
02305 }

const TGraph * NuSystFitter::TwoSidedNeutrinoContourFinder ( const NuMMParameters inputPars  )  [virtual]

Definition at line 276 of file NuSystFitter.cxx.

00277 {
00278   fsn2PointsForNeutrinoContour = fsn2PointsForNeutrinoContourPhysical;
00279   std::vector<std::pair<Double_t,Double_t> > contPointsPhys =
00280     this->NeutrinoContourFinder(inputPars);
00281 
00282   fsn2PointsForNeutrinoContour = fsn2PointsForNeutrinoContourUnphysical;
00283   std::vector<std::pair<Double_t,Double_t> > contPointsUnphys =
00284     this->NeutrinoContourFinder(inputPars);
00285 
00286   //Make a single vector of points
00287   std::vector<pair<Double_t,Double_t> > contourPoints;
00288   for (vector<pair<Double_t, Double_t> >::const_iterator it =
00289          contPointsPhys.begin();
00290        it != contPointsPhys.end();
00291        ++it){
00292     pair<Double_t,Double_t> currentPoint((*it).first,(*it).second);
00293     contourPoints.push_back(currentPoint);
00294   }
00295   for (vector<pair<Double_t, Double_t> >::reverse_iterator
00296          itr = contPointsUnphys.rbegin();
00297        itr != contPointsUnphys.rend();
00298        ++itr){
00299     pair<Double_t,Double_t> currentPoint((*itr).first,(*itr).second);
00300     contourPoints.push_back(currentPoint);
00301   }
00302 
00303 
00304   //Put the contour points in a TGraph
00305   TGraph* contour = new TGraph(contourPoints.size());
00306   Int_t point=0;
00307   for (vector<pair<Double_t, Double_t> >::const_iterator it =
00308          contourPoints.begin();
00309        it != contourPoints.end();
00310        ++it){
00311     contour->SetPoint(point,(*it).first,(*it).second);
00312     cout << "Set contour point " << point
00313          << ", " << (*it).first
00314          << ", " << (*it).second
00315          << endl;
00316     ++point;
00317   }
00318   return contour;
00319 }

virtual double NuSystFitter::Up (  )  const [inline, virtual]

Definition at line 24 of file NuSystFitter.h.

References fUp.

00024 {return fUp;}//Chi2 errors; 0.5 for likelihood


Member Data Documentation

Int_t NuSystFitter::fNumContourPoints [protected]

Definition at line 281 of file NuSystFitter.h.

Referenced by NumberOfContourPoints(), and PlotContours().

Int_t NuSystFitter::fQuietMode [protected]

Definition at line 279 of file NuSystFitter.h.

Referenced by BatchModeOn(), Chi2SliceDm2(), Chi2SliceDm2bar(), Chi2SliceSin2bar(), Chi2ValleyBar(), Chi2ValleyBarDm2Bar(), Chi2ValleyDm2(), Chi2ValleySn2(), FillLikelihoodSurfaceCPT(), FillLikelihoodSurfaceGeneric(), FillLikelihoodSurfaceNSIDm2Eps(), FillLikelihoodSurfaceNSIDm2Sn2(), FillLikelihoodSurfaceNSISn2Eps(), FillLikelihoodSurfaceTransDm2BarGMu(), FillLikelihoodSurfaceTransDm2BarSn2Bar(), FillLikelihoodSurfaceTransDm2GMu(), FillLikelihoodSurfaceTransDm2GMuBoth(), FillLikelihoodSurfaceTransDm2Sn2(), FillLikelihoodSurfaceTransDm2Sn2Both(), FillLikelihoodSurfaceTransSn2BarGMu(), FillLikelihoodSurfaceTransSn2GMu(), FillLikelihoodSurfaceTransSn2GMuBoth(), Likelihood(), LikelihoodSurfaceBar(), Minimise(), QuietModeOff(), and QuietModeOn().

std::vector<NuMMRun*>* NuSystFitter::fRuns [protected]

Definition at line 285 of file NuSystFitter.h.

Referenced by NuABFitter::Archive(), Likelihood(), push_back(), QuietModeOff(), and QuietModeOn().

std::vector<Double_t> NuSystFitter::fsn2PointsForNeutrinoContour [protected]

Definition at line 282 of file NuSystFitter.h.

Referenced by Sn2PointsForNeutrinoContour().

std::vector<Double_t> NuSystFitter::fsn2PointsForNeutrinoContourPhysical [protected]

Definition at line 283 of file NuSystFitter.h.

std::vector<Double_t> NuSystFitter::fsn2PointsForNeutrinoContourUnphysical [protected]

Definition at line 284 of file NuSystFitter.h.

Double_t NuSystFitter::fUp [protected]

Definition at line 280 of file NuSystFitter.h.

Referenced by SetErrorDef(), and Up().


The documentation for this class was generated from the following files:
Generated on Tue Mar 10 23:01:59 2015 for loon by  doxygen 1.4.7