NuSystFitter Class Reference

#include <NuSystFitter.h>

Inheritance diagram for NuSystFitter:

NuABFitter List of all members.

Public Member Functions

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

Protected Member Functions

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

Protected Attributes

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

Detailed Description

Definition at line 19 of file NuSystFitter.h.


Constructor & Destructor Documentation

NuSystFitter::NuSystFitter (  ) 

Definition at line 34 of file NuSystFitter.cxx.

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

NuSystFitter::~NuSystFitter (  )  [virtual]

Definition at line 61 of file NuSystFitter.cxx.

00062 {
00063 }


Member Function Documentation

void NuSystFitter::BatchModeOn (  )  [virtual]

Supresses much of the usual output.

Definition at line 679 of file NuSystFitter.cxx.

References fQuietMode, and QuietModeOn().

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

00680 {
00681   QuietModeOn();
00682   fQuietMode = 2;
00683 }

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

Definition at line 2608 of file NuSystFitter.cxx.

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

Referenced by FCMinimise().

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

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

Definition at line 2933 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNSI().

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

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

Definition at line 2768 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNu().

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

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

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

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

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

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

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

Definition at line 2734 of file NuSystFitter.cxx.

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

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

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

Definition at line 3388 of file NuSystFitter.cxx.

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

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

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

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

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

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

Referenced by GridSearchNu().

01389 {
01390   mmPars.ContourForNuMus();
01391   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01392 }

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

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

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

01398 {
01399   mmPars.ContourForNuBars();
01400   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01401 }

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

Definition at line 1407 of file NuSystFitter.cxx.

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

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

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

Definition at line 2524 of file NuSystFitter.cxx.

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

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

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

Definition at line 2494 of file NuSystFitter.cxx.

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

02497 {
02498   mmPars.ContourForDecay(ContourPair);
02499   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02500 }

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

Definition at line 2504 of file NuSystFitter.cxx.

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

02507 {
02508   mmPars.ContourForDecoherence(ContourPair);
02509   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02510 }

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

Definition at line 2514 of file NuSystFitter.cxx.

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

02517 {
02518   mmPars.ContourForDecoherenceBar(ContourPair);
02519   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02520 }

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

Definition at line 2311 of file NuSystFitter.cxx.

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

02313 {
02314   mmPars.ContourForDm2BarSn2BarGMu3D(); 
02315   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02316 }

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

Definition at line 1442 of file NuSystFitter.cxx.

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

01444 {
01445   mmPars.ContourForDm2GMuCMu3D();  
01446   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01447 }

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

Definition at line 2320 of file NuSystFitter.cxx.

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

02322 {
02323   mmPars.ContourForDm2Sn2GMu3D(); 
02324   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02325 }

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

Definition at line 2329 of file NuSystFitter.cxx.

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

02331 {
02332   mmPars.ContourForDm2Sn2GMu3DBoth(); 
02333   return FillLikelihoodSurfaceGeneric(surface, mmPars);
02334 }

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

Definition at line 1322 of file NuSystFitter.cxx.

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

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

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

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

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

Definition at line 1434 of file NuSystFitter.cxx.

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

01436 {                                               
01437   mmPars.ContourForGMuCMuMu();
01438   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01439 }

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

Definition at line 1416 of file NuSystFitter.cxx.

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

01418 {                                               
01419   mmPars.ContourForGMuGTau();
01420   return FillLikelihoodSurfaceGeneric(surface, mmPars);
01421 }

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

Definition at line 1425 of file NuSystFitter.cxx.

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

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

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

Definition at line 2567 of file NuSystFitter.cxx.

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

02567                                                                                    {
02568   mmPars.ContourForLED();
02569   return FillLikelihoodSurfaceGeneric(surface, mmPars);    
02570 }

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

Definition at line 2338 of file NuSystFitter.cxx.

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

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

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

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

Definition at line 2358 of file NuSystFitter.cxx.

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

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

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

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

Definition at line 2348 of file NuSystFitter.cxx.

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

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

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

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

Definition at line 2373 of file NuSystFitter.cxx.

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

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

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

Definition at line 2453 of file NuSystFitter.cxx.

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

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

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

Definition at line 2412 of file NuSystFitter.cxx.

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

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

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

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

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

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

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

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

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

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

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

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

Definition at line 1451 of file NuSystFitter.cxx.

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

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

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

Definition at line 1688 of file NuSystFitter.cxx.

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

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

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

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

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

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

Definition at line 1650 of file NuSystFitter.cxx.

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

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

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

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

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

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

Definition at line 2204 of file NuSystFitter.cxx.

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

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

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

Definition at line 2103 of file NuSystFitter.cxx.

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

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

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

Definition at line 1772 of file NuSystFitter.cxx.

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

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

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

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

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

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

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

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

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

Definition at line 2171 of file NuSystFitter.cxx.

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

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

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

Definition at line 2070 of file NuSystFitter.cxx.

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

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

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

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

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

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

Definition at line 1530 of file NuSystFitter.cxx.

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

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

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

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

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

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

Definition at line 1610 of file NuSystFitter.cxx.

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

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

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

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

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

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

Definition at line 2237 of file NuSystFitter.cxx.

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

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

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

Definition at line 2136 of file NuSystFitter.cxx.

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

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

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

Definition at line 2641 of file NuSystFitter.cxx.

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

Referenced by FCMinimise().

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

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

Definition at line 2961 of file NuSystFitter.cxx.

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

Referenced by FCMinimiseNSI().

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

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

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

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

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

Definition at line 2694 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearch(), and FineGridSearch().

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

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

Definition at line 3019 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearchNSI(), and FineGridSearchNSI().

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

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

Definition at line 2854 of file NuSystFitter.cxx.

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

Referenced by CoarseGridSearchNu(), and FineGridSearchNu().

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

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         cout << "stats chi2 = " << likelihood << endl;
00669     likelihood += mmPars.PenaltyTerm();
00670         cout << "stats + penalty = " << likelihood << endl;
00671     if (!fQuietMode){
00672 //        cout << " + penalty = " << likelihood << endl;
00673     }
00674 
00675     return likelihood;
00676 }

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

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

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

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 Chi2ValleyBar(), Chi2ValleyBarDm2Bar(), Chi2ValleyDm2(), Chi2ValleySn2(), 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(), LikelihoodSurfaceBar(), NuFCFitterNSINubar::MinuitFit(), NuFCFitterNSINu::MinuitFit(), NuFCFitterNSI::MinuitFit(), NuFCFitter::MinuitFit(), 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(), 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 697 of file NuSystFitter.cxx.

References fQuietMode, and fRuns.

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

void NuSystFitter::QuietModeOn (  )  [virtual]

Definition at line 686 of file NuSystFitter.cxx.

References fQuietMode, and fRuns.

Referenced by BatchModeOn().

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

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

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

Referenced by FillCPT().

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

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

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

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

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

Definition at line 276 of file NuSystFitter.cxx.

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

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

Definition at line 24 of file NuSystFitter.h.

References fUp.

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


Member Data Documentation

Int_t NuSystFitter::fNumContourPoints [protected]

Definition at line 281 of file NuSystFitter.h.

Referenced by NumberOfContourPoints(), and PlotContours().

Int_t NuSystFitter::fQuietMode [protected]

Definition at line 279 of file NuSystFitter.h.

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

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

Definition at line 285 of file NuSystFitter.h.

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

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

Definition at line 282 of file NuSystFitter.h.

Referenced by Sn2PointsForNeutrinoContour().

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

Definition at line 283 of file NuSystFitter.h.

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

Definition at line 284 of file NuSystFitter.h.

Double_t NuSystFitter::fUp [protected]

Definition at line 280 of file NuSystFitter.h.

Referenced by SetErrorDef(), and Up().


The documentation for this class was generated from the following files:
Generated on Mon Aug 11 01:06:38 2014 for loon by  doxygen 1.4.7