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, 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 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 3197 of file NuSystFitter.cxx.

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

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

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 2611 of file NuSystFitter.cxx.

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

Referenced by FCMinimise().

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

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

Definition at line 2936 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNSI().

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

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

Definition at line 2771 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNu().

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

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 2577 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().

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

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

Definition at line 2896 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().

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

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

Definition at line 2737 of file NuSystFitter.cxx.

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

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

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

Definition at line 3391 of file NuSystFitter.cxx.

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

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

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 3320 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().

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

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,
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 2565 of file NuSystFitter.cxx.

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

02570                                                                                         {
02571   mmPars.ContourForLED(ContourPair);
02572   return FillLikelihoodSurfaceGeneric(surface, mmPars, Ix, Iy, Nx, Ny, mmPars_surfaces, converge_surface, par_index, n_iter, iter);    
02573 }

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 2644 of file NuSystFitter.cxx.

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

Referenced by FCMinimise().

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

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

Definition at line 2964 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNSI().

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

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

Definition at line 2804 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().

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

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

Definition at line 2697 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearch(), and FineGridSearch().

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

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

Definition at line 3022 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearchNSI(), and FineGridSearchNSI().

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

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

Definition at line 2857 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearchNu(), and FineGridSearchNu().

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

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::CalculateChi2Best1D(), 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 3253 of file NuSystFitter.cxx.

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

Referenced by FillCPT().

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

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 285 of file NuSystFitter.h.

Referenced by NumberOfContourPoints(), and PlotContours().

Int_t NuSystFitter::fQuietMode [protected]

Definition at line 283 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 289 of file NuSystFitter.h.

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

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

Definition at line 286 of file NuSystFitter.h.

Referenced by Sn2PointsForNeutrinoContour().

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

Definition at line 287 of file NuSystFitter.h.

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

Definition at line 288 of file NuSystFitter.h.

Double_t NuSystFitter::fUp [protected]

Definition at line 284 of file NuSystFitter.h.

Referenced by SetErrorDef(), and Up().


The documentation for this class was generated from the following files:
Generated on Fri Jul 10 22:45:28 2015 for loon by  doxygen 1.4.7