MaxLikelihoodFit Class Reference

#include <MaxLikelihoodFit.h>

Inheritance diagram for MaxLikelihoodFit:

Public Member Functions
virtual void	applyOptions (const boost::program_options::variables_map &vm)
	MaxLikelihoodFit ()
virtual const std::string &	name () const
virtual void	setNToys (const int)
virtual void	setToyNumber (const int)
	~MaxLikelihoodFit ()
Public Member Functions inherited from FitterAlgoBase
void	applyOptionsBase (const boost::program_options::variables_map &vm)
	FitterAlgoBase (const char *title="<FillMe> specific options")
virtual bool	run (RooWorkspace *w, RooStats::ModelConfig *mc_s, RooStats::ModelConfig *mc_b, RooAbsData &data, double &limit, double &limitErr, const double *hint)
Public Member Functions inherited from LimitAlgo
virtual void	applyDefaultOptions ()
	LimitAlgo ()
	LimitAlgo (const char *desc)
const boost::program_options::options_description &	options () const

Protected Member Functions
void	createFitResultTrees (const RooStats::ModelConfig &)
void	getNormalizations (RooAbsPdf *pdf, const RooArgSet &obs, RooArgSet &out)
virtual bool	runSpecific (RooWorkspace *w, RooStats::ModelConfig *mc_s, RooStats::ModelConfig *mc_b, RooAbsData &data, double &limit, double &limitErr, const double *hint)
void	setFitResultTrees (const RooArgSet *, double *)
Protected Member Functions inherited from FitterAlgoBase
RooFitResult *	doFit (RooAbsPdf &pdf, RooAbsData &data, RooRealVar &r, const RooCmdArg &constrain, bool doHesse=true, int ndim=1, bool reuseNLL=false)
RooFitResult *	doFit (RooAbsPdf &pdf, RooAbsData &data, const RooArgList &rs, const RooCmdArg &constrain, bool doHesse=true, int ndim=1, bool reuseNLL=false)
double	findCrossing (CascadeMinimizer &minim, RooAbsReal &nll, RooRealVar &r, double level, double rStart, double rBound)

Protected Attributes
int	currentToy_
std::auto_ptr< TFile >	fitOut
int	fitStatus_
double *	globalObservables_
double	mu_
double	nll_bonly_
double	nll_nll0_
double	nll_sb_
int	nToys
double *	nuisanceParameters_
int	numbadnll_
TTree *	t_fit_b_
TTree *	t_fit_sb_
Protected Attributes inherited from FitterAlgoBase
std::auto_ptr< RooAbsReal >	nll
Protected Attributes inherited from LimitAlgo
boost::program_options::options_description	options_

Static Protected Attributes
static std::string	backgroundPdfNames_ = "shapeBkg*"
static bool	justFit_ = false
static bool	makePlots_ = false
static std::string	minos_ = "poi"
static std::string	name_ = ""
static bool	noErrors_ = false
static std::string	out_ = "."
static float	rebinFactor_ = 1.0
static bool	reuseParams_ = false
static bool	saveNormalizations_ = false
static std::string	signalPdfNames_ = "shapeSig*"
Static Protected Attributes inherited from FitterAlgoBase
static bool	do95_ = false
static bool	keepFailures_ = false
static int	maxFailedSteps_ = 5
static std::string	minimizerAlgo_ = "Minuit2"
static std::string	minimizerAlgoForMinos_ = "Minuit2,simplex"
static int	minimizerStrategy_ = 1
static int	minimizerStrategyForMinos_ = 0
static float	minimizerTolerance_ = 1e-2
static float	minimizerToleranceForMinos_ = 1e-4
static float	nllValue_ = std::numeric_limits<float>::quiet_NaN()
static float	preFitValue_ = 1.0
static bool	robustFit_ = false
static bool	saveNLL_ = false
static float	stepSize_ = 0.1

Detailed Description

Do a ML fit of the data with background and signal+background hypothesis and print out diagnostics plots

Author

Giovanni Petrucciani (UCSD)

Definition at line 13 of file MaxLikelihoodFit.h.

Constructor & Destructor Documentation

MaxLikelihoodFit::MaxLikelihoodFit

(

)

Definition at line 43 of file MaxLikelihoodFit.cc.

References backgroundPdfNames_, fitStatus_, minos_, mu_, nll_bonly_, nll_nll0_, nll_sb_, nToys, numbadnll_, LimitAlgo::options_, out_, rebinFactor_, and signalPdfNames_.

                                   :
    FitterAlgoBase("MaxLikelihoodFit specific options")
{
    options_.add_options()
        ("minos",              boost::program_options::value<std::string>(&minos_)->default_value(minos_), "Compute MINOS errors for: 'none', 'poi', 'all'")
        ("out",                boost::program_options::value<std::string>(&out_)->default_value(out_), "Directory to put output in")
        ("plots",              "Make plots")
        ("rebinFactor",        boost::program_options::value<float>(&rebinFactor_)->default_value(rebinFactor_), "Rebin by this factor before plotting (does not affect fitting!)")
        ("signalPdfNames",     boost::program_options::value<std::string>(&signalPdfNames_)->default_value(signalPdfNames_), "Names of signal pdfs in plots (separated by ,)")
        ("backgroundPdfNames", boost::program_options::value<std::string>(&backgroundPdfNames_)->default_value(backgroundPdfNames_), "Names of background pdfs in plots (separated by ',')")
        ("saveNormalizations",  "Save post-fit normalizations of all components of the pdfs")
        ("justFit",  "Just do the S+B fit, don't do the B-only one, don't save output file")
        ("noErrors",  "Don't compute uncertainties on the best fit value")
        ("initFromBonly",  "Use the vlaues of the nuisance parameters from the background only fit as the starting point for the s+b fit")
   ;

    // setup a few defaults
    nToys=0; fitStatus_=0; mu_=0;numbadnll_=-1;nll_nll0_=-1; nll_bonly_=-1;nll_sb_=-1;
}

MaxLikelihoodFit::~MaxLikelihoodFit

(

)

Definition at line 63 of file MaxLikelihoodFit.cc.

References globalObservables_, and nuisanceParameters_.

                                   {
   // delete the Arrays used to fill the trees;
   delete globalObservables_;
   delete nuisanceParameters_;
}

Member Function Documentation

void MaxLikelihoodFit::applyOptions ( const boost::program_options::variables_map &  vm )

virtual

Reimplemented from LimitAlgo.

Definition at line 75 of file MaxLikelihoodFit.cc.

References FitterAlgoBase::applyOptionsBase(), justFit_, makePlots_, name_, noErrors_, out_, reuseParams_, and saveNormalizations_.

{
    applyOptionsBase(vm);
    makePlots_ = vm.count("plots");
    name_ = vm["name"].defaulted() ?  std::string() : vm["name"].as<std::string>();
    saveNormalizations_  = vm.count("saveNormalizations");
    justFit_  = vm.count("justFit");
    noErrors_ = vm.count("noErrors");
    reuseParams_ = vm.count("initFromBonly");
    if (justFit_) { out_ = "none"; makePlots_ = false; saveNormalizations_ = false; reuseParams_ = false;}
    // For now default this to true;
}

void MaxLikelihoodFit::createFitResultTrees ( const RooStats::ModelConfig &  mc )

protected

Definition at line 402 of file MaxLikelihoodFit.cc.

References a, prof2calltree::count, gather_cfg::cout, fitStatus_, globalObservables_, mu_, name(), nll_bonly_, nll_nll0_, nll_sb_, nuisanceParameters_, numbadnll_, t_fit_b_, and t_fit_sb_.

Referenced by runSpecific().

                                                                        {

     // Initiate the arrays to store parameters

     // create TTrees to store fit results:
     t_fit_b_  = new TTree("tree_fit_b","tree_fit_b");
     t_fit_sb_ = new TTree("tree_fit_sb","tree_fit_sb");

         t_fit_b_->Branch("fit_status",&fitStatus_,"fit_status/Int_t");
     t_fit_sb_->Branch("fit_status",&fitStatus_,"fit_status/Int_t");

     t_fit_b_->Branch("mu",&mu_,"mu/Double_t");
     t_fit_sb_->Branch("mu",&mu_,"mu/Double_t");

     t_fit_b_->Branch("numbadnll",&numbadnll_,"numbadnll/Int_t");
     t_fit_sb_->Branch("numbadnll",&numbadnll_,"numbadnll/Int_t");

     t_fit_b_->Branch("nll_min",&nll_bonly_,"nll_min/Double_t");
     t_fit_sb_->Branch("nll_min",&nll_sb_,"nll_min/Double_t");

     t_fit_sb_->Branch("nll_nll0",&nll_nll0_,"nll_nll0/Double_t");

         // fill the maps for the nuisances, and global observables
         const RooArgSet *cons = mc.GetGlobalObservables();
         const RooArgSet *nuis = mc.GetNuisanceParameters();

     globalObservables_ = new double[cons->getSize()];
     nuisanceParameters_= new double[nuis->getSize()];
        
     int count=0; 
         TIterator* iter_c(cons->createIterator());
         for (TObject *a = iter_c->Next(); a != 0; a = iter_c->Next()) { 
                 RooRealVar *rrv = dynamic_cast<RooRealVar *>(a);        
         std::string name = rrv->GetName();
         globalObservables_[count]=0;
         t_fit_sb_->Branch(name.c_str(),&(globalObservables_[count]),Form("%s/Double_t",name.c_str()));
         t_fit_b_->Branch(name.c_str(),&(globalObservables_[count]),Form("%s/Double_t",name.c_str()));
         count++;
         }
         
     count = 0;
         TIterator* iter_n(nuis->createIterator());
         for (TObject *a = iter_n->Next(); a != 0; a = iter_n->Next()) { 
                 RooRealVar *rrv = dynamic_cast<RooRealVar *>(a);        
         std::string name = rrv->GetName();
         nuisanceParameters_[count] = 0;
         t_fit_sb_->Branch(name.c_str(),&(nuisanceParameters_[count])),Form("%s/Double_t",name.c_str());
         t_fit_b_->Branch(name.c_str(),&(nuisanceParameters_[count]),Form("%s/Double_t",name.c_str()));
         count++;
         }
    std::cout << "Created Branches" <<std::endl;
         return;    
}

void MaxLikelihoodFit::getNormalizations ( RooAbsPdf *  pdf, const RooArgSet &  obs, RooArgSet &  out  )

protected

Definition at line 355 of file MaxLikelihoodFit.cc.

References Clusterizer1DCommons::add(), i, list(), n, and parseEventContent::prod.

Referenced by runSpecific().

                                                                                             {
    RooSimultaneous *sim = dynamic_cast<RooSimultaneous *>(pdf);
    if (sim != 0) {
        RooAbsCategoryLValue &cat = const_cast<RooAbsCategoryLValue &>(sim->indexCat());
        for (int i = 0, n = cat.numBins((const char *)0); i < n; ++i) {
            cat.setBin(i);
            RooAbsPdf *pdfi = sim->getPdf(cat.getLabel());
            if (pdfi) getNormalizations(pdfi, obs, out);
        }        
        return;
    }
    RooProdPdf *prod = dynamic_cast<RooProdPdf *>(pdf);
    if (prod != 0) {
        RooArgList list(prod->pdfList());
        for (int i = 0, n = list.getSize(); i < n; ++i) {
            RooAbsPdf *pdfi = (RooAbsPdf *) list.at(i);
            if (pdfi->dependsOn(obs)) getNormalizations(pdfi, obs, out);
        }
        return;
    }
    RooAddPdf *add = dynamic_cast<RooAddPdf *>(pdf);
    if (add != 0) {
        RooArgList list(add->coefList());
        for (int i = 0, n = list.getSize(); i < n; ++i) {
            RooAbsReal *coeff = (RooAbsReal *) list.at(i);
            out.addOwned(*(new RooConstVar(coeff->GetName(), "", coeff->getVal())));
        }
        return;
    }
}

virtual const std::string& MaxLikelihoodFit::name ( void  ) const

inlinevirtual

Implements LimitAlgo.

Definition at line 16 of file MaxLikelihoodFit.h.

Referenced by BeautifulSoup.Tag::_invert(), createFitResultTrees(), and setFitResultTrees().

                                         {
    static const std::string name("MaxLikelihoodFit");
    return name;
  }

bool MaxLikelihoodFit::runSpecific ( RooWorkspace *  w, RooStats::ModelConfig *  mc_s, RooStats::ModelConfig *  mc_b, RooAbsData &  data, double &  limit, double &  limitErr, const double *  hint  )

protectedvirtual

Implements FitterAlgoBase.

Definition at line 88 of file MaxLikelihoodFit.cc.

References backgroundPdfNames_, alignmentValidation::c1, Combine::commitPoint(), corr, gather_cfg::cout, createFitResultTrees(), currentToy_, data, FitterAlgoBase::do95_, FitterAlgoBase::doFit(), fitOut, fitStatus_, getNormalizations(), globalObservables_, justFit_, utils::makeNuisancePdf(), utils::makePlots(), makePlots_, FitterAlgoBase::minimizerStrategy_, minos_, mu_, name_, FitterAlgoBase::nll, nll_bonly_, nll_nll0_, nll_sb_, noErrors_, nToys, nuisanceParameters_, numbadnll_, out_, RecoTauValidation_cfi::plots, FitterAlgoBase::preFitValue_, alignCSCRings::r, rebinFactor_, reuseParams_, saveNormalizations_, setFitResultTrees(), signalPdfNames_, t_fit_b_, t_fit_sb_, utils::tdrStyle(), and withSystematics.

                                                                                                                                                                               {

  if (reuseParams_ && minos_!="none"){
    std::cout << "Cannot reuse b-only fit params when running minos. Parameters will be reset when running S+B fit"<<std::endl;
    reuseParams_=false;
  }

  if (!justFit_ && out_ != "none"){
    if (currentToy_ < 1){
        fitOut.reset(TFile::Open((out_+"/mlfit"+name_+".root").c_str(), "RECREATE")); 
        createFitResultTrees(*mc_s);
    }
  }

  RooRealVar *r = dynamic_cast<RooRealVar *>(mc_s->GetParametersOfInterest()->first());

  TCanvas *c1 = 0;
  if (makePlots_) {
      utils::tdrStyle();
      c1 = new TCanvas("c1","c1");
  }

  // Make pre-plots before the fit
  r->setVal(preFitValue_);
  if (makePlots_) {
      std::vector<RooPlot *> plots = utils::makePlots(*mc_s->GetPdf(), data, signalPdfNames_.c_str(), backgroundPdfNames_.c_str(), rebinFactor_);
      for (std::vector<RooPlot *>::iterator it = plots.begin(), ed = plots.end(); it != ed; ++it) {
          (*it)->Draw(); 
          c1->Print((out_+"/"+(*it)->GetName()+"_prefit.png").c_str());
          if (fitOut.get() && currentToy_< 1) fitOut->WriteTObject(*it, (std::string((*it)->GetName())+"_prefit").c_str());
      }
  }


  // Determine pre-fit values of nuisance parameters
  if (currentToy_ < 1){
    const RooArgSet *nuis      = mc_s->GetNuisanceParameters();
    const RooArgSet *globalObs = mc_s->GetGlobalObservables();
    if (!justFit_ && nuis && globalObs ) {
      std::auto_ptr<RooAbsPdf> nuisancePdf(utils::makeNuisancePdf(*mc_s));
      std::auto_ptr<RooDataSet> globalData(new RooDataSet("globalData","globalData", *globalObs));
      globalData->add(*globalObs);
      RooFitResult *res_prefit = 0;
      {     
            CloseCoutSentry sentry(verbose < 2);
            res_prefit = nuisancePdf->fitTo(*globalData,
            RooFit::Save(1),
            RooFit::Minimizer(ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str(), ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str()),
            RooFit::Strategy(minimizerStrategy_),
            RooFit::Minos(minos_ == "all")
            );
      }
      if (fitOut.get() ) fitOut->WriteTObject(res_prefit, "nuisances_prefit_res");
      if (fitOut.get() ) fitOut->WriteTObject(nuis->snapshot(), "nuisances_prefit");

      nuisancePdf.reset();
      globalData.reset();
      delete res_prefit;

    } else if (nuis) {
      if (fitOut.get() ) fitOut->WriteTObject(nuis->snapshot(), "nuisances_prefit");
    }
  }
  
  RooFitResult *res_b = 0, *res_s = 0;
  const RooCmdArg &constCmdArg_s = withSystematics  ? RooFit::Constrain(*mc_s->GetNuisanceParameters()) : RooFit::NumCPU(1); // use something dummy 
  const RooCmdArg &minosCmdArg = minos_ == "poi" ?  RooFit::Minos(*mc_s->GetParametersOfInterest())   : RooFit::Minos(minos_ != "none"); 
  w->loadSnapshot("clean");
  r->setVal(0.0); r->setConstant(true);

  // Setup Nll before calling fits;
  if (currentToy_<1) nll.reset(mc_s->GetPdf()->createNLL(data,constCmdArg_s,RooFit::Extended(mc_s->GetPdf()->canBeExtended())));
  // Get the nll value on the prefit
  double nll0 = nll->getVal();

  if (justFit_) { 
    // skip b-only fit
  } else if (minos_ != "all") {
    RooArgList minos; 
    res_b = doFit(*mc_s->GetPdf(), data, minos, constCmdArg_s, /*hesse=*/true,/*reuseNLL*/ true); 
    nll_bonly_=nll->getVal()-nll0;   
  } else {
    CloseCoutSentry sentry(verbose < 2);
    res_b = mc_s->GetPdf()->fitTo(data, 
            RooFit::Save(1), 
            RooFit::Minimizer(ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str(), ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str()), 
            RooFit::Strategy(minimizerStrategy_),
            RooFit::Extended(mc_s->GetPdf()->canBeExtended()), 
            constCmdArg_s, minosCmdArg
            );
    if (res_b) nll_bonly_ = nll->getVal() - nll0;

  }

  if (res_b) { 
      if (verbose > 1) res_b->Print("V");
      if (fitOut.get()) {
     if (currentToy_< 1)    fitOut->WriteTObject(res_b,"fit_b");
     setFitResultTrees(mc_s->GetNuisanceParameters(),nuisanceParameters_);
     setFitResultTrees(mc_s->GetGlobalObservables(),globalObservables_);
     fitStatus_ = res_b->status();
      }
      numbadnll_=res_b->numInvalidNLL();

      if (makePlots_) {
          std::vector<RooPlot *> plots = utils::makePlots(*mc_b->GetPdf(), data, signalPdfNames_.c_str(), backgroundPdfNames_.c_str(), rebinFactor_);
          for (std::vector<RooPlot *>::iterator it = plots.begin(), ed = plots.end(); it != ed; ++it) {
              c1->cd(); (*it)->Draw(); 
              c1->Print((out_+"/"+(*it)->GetName()+"_fit_b.png").c_str());
              if (fitOut.get() && currentToy_< 1) fitOut->WriteTObject(*it, (std::string((*it)->GetName())+"_fit_b").c_str());
          }
      }

      if (saveNormalizations_ && currentToy_<1) {
          RooArgSet *norms = new RooArgSet();
          norms->setName("norm_fit_b");
          getNormalizations(mc_s->GetPdf(), *mc_s->GetObservables(), *norms);
          if (fitOut.get()) fitOut->WriteTObject(norms, "norm_fit_b");
      delete norms;
      }

      if (makePlots_ && currentToy_<1)  {
      TH2 *corr = res_b->correlationHist();
          c1->SetLeftMargin(0.25);  c1->SetBottomMargin(0.25);
          corr->SetTitle("Correlation matrix of fit parameters");
          gStyle->SetPaintTextFormat(res_b->floatParsFinal().getSize() > 10 ? ".1f" : ".2f");
          gStyle->SetOptStat(0);
          corr->SetMarkerSize(res_b->floatParsFinal().getSize() > 10 ? 2 : 1);
          corr->Draw("COLZ TEXT");
          c1->Print((out_+"/covariance_fit_b.png").c_str());
          c1->SetLeftMargin(0.16);  c1->SetBottomMargin(0.13);
          if (fitOut.get()) fitOut->WriteTObject(corr, "covariance_fit_b");
      }
  }
  else {
    fitStatus_=-1;
    numbadnll_=-1;  
  }
  mu_=r->getVal();
  if (t_fit_b_) t_fit_b_->Fill();
  // no longer need res_b
  delete res_b;

  if (!reuseParams_) w->loadSnapshot("clean"); // Reset, also ensures nll_prefit is same in call to doFit for b and s+b
  r->setVal(preFitValue_); r->setConstant(false); 
  if (minos_ != "all") {
    RooArgList minos; if (minos_ == "poi") minos.add(*r);
    res_s = doFit(*mc_s->GetPdf(), data, minos, constCmdArg_s, /*hesse=*/!noErrors_,/*reuseNLL*/ true); 
    nll_sb_ = nll->getVal()-nll0;
  } else {
    CloseCoutSentry sentry(verbose < 2);
    res_s = mc_s->GetPdf()->fitTo(data, 
            RooFit::Save(1), 
            RooFit::Minimizer(ROOT::Math::MinimizerOptions::DefaultMinimizerType().c_str(), ROOT::Math::MinimizerOptions::DefaultMinimizerAlgo().c_str()), 
            RooFit::Strategy(minimizerStrategy_),
            RooFit::Extended(mc_s->GetPdf()->canBeExtended()), 
            constCmdArg_s, minosCmdArg
            );
    if (res_s) nll_sb_= nll->getVal()-nll0;

  }
  if (res_s) { 
      limit    = r->getVal();
      limitErr = r->getError();
      if (verbose > 1) res_s->Print("V");
      if (fitOut.get()){
     if (currentToy_<1) fitOut->WriteTObject(res_s, "fit_s");

     setFitResultTrees(mc_s->GetNuisanceParameters(),nuisanceParameters_);
     setFitResultTrees(mc_s->GetGlobalObservables(),globalObservables_);
     fitStatus_ = res_s->status();
         numbadnll_ = res_s->numInvalidNLL();

     // Additionally store the nll_sb - nll_bonly (=0.5*q0)
     nll_nll0_ =  nll_sb_ -  nll_bonly_;
      }

      if (makePlots_) {
          std::vector<RooPlot *> plots = utils::makePlots(*mc_s->GetPdf(), data, signalPdfNames_.c_str(), backgroundPdfNames_.c_str(), rebinFactor_);
          for (std::vector<RooPlot *>::iterator it = plots.begin(), ed = plots.end(); it != ed; ++it) {
              c1->cd(); (*it)->Draw(); 
              c1->Print((out_+"/"+(*it)->GetName()+"_fit_s.png").c_str());
              if (fitOut.get() && currentToy_< 1) fitOut->WriteTObject(*it, (std::string((*it)->GetName())+"_fit_s").c_str());
          }
      }

      if (saveNormalizations_&& currentToy_< 1) {
          RooArgSet *norms = new RooArgSet();
          norms->setName("norm_fit_s");
          getNormalizations(mc_s->GetPdf(), *mc_s->GetObservables(), *norms);
          if (fitOut.get() ) fitOut->WriteTObject(norms, "norm_fit_s");
      delete norms;
      }

      if (makePlots_&& currentToy_< 1)  {
          TH2 *corr = res_s->correlationHist();
          c1->SetLeftMargin(0.25);  c1->SetBottomMargin(0.25);
          corr->SetTitle("Correlation matrix of fit parameters");
          gStyle->SetPaintTextFormat(res_s->floatParsFinal().getSize() > 10 ? ".1f" : ".2f");
          gStyle->SetOptStat(0);
          corr->SetMarkerSize(res_s->floatParsFinal().getSize() > 10 ? 2 : 1);
          corr->Draw("COLZ TEXT");
          c1->Print((out_+"/covariance_fit_s.png").c_str());
          c1->SetLeftMargin(0.16);  c1->SetBottomMargin(0.13);
          if (fitOut.get() ) fitOut->WriteTObject(corr, "covariance_fit_s");
      }
  }  else {
    fitStatus_=-1;
    numbadnll_=-1;
    nll_nll0_ = -1;
  }
  mu_=r->getVal();
  if (t_fit_sb_) t_fit_sb_->Fill();

  if (res_s) {
      RooRealVar *rf = dynamic_cast<RooRealVar*>(res_s->floatParsFinal().find(r->GetName()));
      double bestFitVal = rf->getVal();

      double hiErr = +(rf->hasRange("err68") ? rf->getMax("err68") - bestFitVal : rf->getAsymErrorHi());
      double loErr = -(rf->hasRange("err68") ? rf->getMin("err68") - bestFitVal : rf->getAsymErrorLo());
      double maxError = std::max<double>(std::max<double>(hiErr, loErr), rf->getError());

      if (fabs(hiErr) < 0.001*maxError) hiErr = -bestFitVal + rf->getMax();
      if (fabs(loErr) < 0.001*maxError) loErr = +bestFitVal - rf->getMin();

      double hiErr95 = +(do95_ && rf->hasRange("err95") ? rf->getMax("err95") - bestFitVal : 0);
      double loErr95 = -(do95_ && rf->hasRange("err95") ? rf->getMin("err95") - bestFitVal : 0);

      limit = bestFitVal;  limitErr = 0;
      if (!noErrors_) Combine::commitPoint(/*expected=*/true, /*quantile=*/0.5);
      limit = bestFitVal - loErr; limitErr = 0;
      if (!noErrors_) Combine::commitPoint(/*expected=*/true, /*quantile=*/0.16);
      limit = bestFitVal + hiErr; limitErr = 0;
      if (!noErrors_) Combine::commitPoint(/*expected=*/true, /*quantile=*/0.84);
      if (do95_ && rf->hasRange("err95") && !noErrors_) {
        limit = rf->getMax("err95"); Combine::commitPoint(/*expected=*/true, /*quantile=*/0.975);
        limit = rf->getMin("err95"); Combine::commitPoint(/*expected=*/true, /*quantile=*/0.025);
      }

      limit = bestFitVal;
      limitErr = maxError;
      std::cout << "\n --- MaxLikelihoodFit ---" << std::endl;
      std::cout << "Best fit " << r->GetName() << ": " << rf->getVal() << "  "<<  -loErr << "/+" << +hiErr << "  (68% CL)" << std::endl;
      if (do95_) {
        std::cout << "         " << r->GetName() << ": " << rf->getVal() << "  "<<  -loErr95 << "/+" << +hiErr95 << "  (95% CL)" << std::endl;
      }
  } else {
      std::cout << "\n --- MaxLikelihoodFit ---" << std::endl;
      std::cout << "Fit failed."  << std::endl;
  }

  if (currentToy_==nToys-1 || nToys==0 ) {
        
        if (fitOut.get()) { 
        fitOut->cd();
        t_fit_sb_->Write(); t_fit_b_->Write();
        fitOut.release()->Close();
    }

  } 
  bool fitreturn = (res_s!=0);
  delete res_s;

  std::cout << "nll S+B -> "<<nll_sb_ << "  nll B -> " << nll_bonly_ <<std::endl;
  return fitreturn;
}

void MaxLikelihoodFit::setFitResultTrees ( const RooArgSet *  args, double *  vals  )

protected

Definition at line 387 of file MaxLikelihoodFit.cc.

References a, prof2calltree::count, and name().

Referenced by runSpecific().

                                                                            {
    
         TIterator* iter(args->createIterator());
     int count=0;
     
         for (TObject *a = iter->Next(); a != 0; a = iter->Next()) { 
                 RooRealVar *rrv = dynamic_cast<RooRealVar *>(a);        
         std::string name = rrv->GetName();
         vals[count]=rrv->getVal();
         count++;
         }
     delete iter;
     return;
}

void MaxLikelihoodFit::setNToys ( const int  iToy )

virtual

Reimplemented from LimitAlgo.

Definition at line 72 of file MaxLikelihoodFit.cc.

References nToys.

                                             {
    nToys = iToy;
}

void MaxLikelihoodFit::setToyNumber ( const int  iToy )

virtual

Reimplemented from LimitAlgo.

Definition at line 69 of file MaxLikelihoodFit.cc.

References currentToy_.

                                                 {
    currentToy_ = iToy;
}

Member Data Documentation

std::string MaxLikelihoodFit::backgroundPdfNames_ = "shapeBkg*"

staticprotected

Definition at line 36 of file MaxLikelihoodFit.h.

Referenced by MaxLikelihoodFit(), and runSpecific().

int MaxLikelihoodFit::currentToy_

protected

Definition at line 39 of file MaxLikelihoodFit.h.

Referenced by runSpecific(), and setToyNumber().

std::auto_ptr<TFile> MaxLikelihoodFit::fitOut

protected

Definition at line 42 of file MaxLikelihoodFit.h.

Referenced by runSpecific().

int MaxLikelihoodFit::fitStatus_

protected

Definition at line 40 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

double* MaxLikelihoodFit::globalObservables_

protected

Definition at line 43 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), runSpecific(), and ~MaxLikelihoodFit().

bool MaxLikelihoodFit::justFit_ = false

staticprotected

Definition at line 32 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

bool MaxLikelihoodFit::makePlots_ = false

staticprotected

Definition at line 34 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

std::string MaxLikelihoodFit::minos_ = "poi"

staticprotected

Definition at line 30 of file MaxLikelihoodFit.h.

Referenced by MaxLikelihoodFit(), and runSpecific().

double MaxLikelihoodFit::mu_

protected

Definition at line 41 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

std::string MaxLikelihoodFit::name_ = ""

staticprotected

Definition at line 28 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

double MaxLikelihoodFit::nll_bonly_

protected

Definition at line 41 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

double MaxLikelihoodFit::nll_nll0_

protected

Definition at line 41 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

double MaxLikelihoodFit::nll_sb_

protected

Definition at line 41 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

bool MaxLikelihoodFit::noErrors_ = false

staticprotected

Definition at line 32 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

int MaxLikelihoodFit::nToys

protected

Definition at line 39 of file MaxLikelihoodFit.h.

Referenced by MaxLikelihoodFit(), runSpecific(), and setNToys().

double* MaxLikelihoodFit::nuisanceParameters_

protected

Definition at line 44 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), runSpecific(), and ~MaxLikelihoodFit().

int MaxLikelihoodFit::numbadnll_

protected

Definition at line 40 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), MaxLikelihoodFit(), and runSpecific().

std::string MaxLikelihoodFit::out_ = "."

staticprotected

Definition at line 33 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), MaxLikelihoodFit(), and runSpecific().

float MaxLikelihoodFit::rebinFactor_ = 1.0

staticprotected

Definition at line 35 of file MaxLikelihoodFit.h.

Referenced by MaxLikelihoodFit(), and runSpecific().

bool MaxLikelihoodFit::reuseParams_ = false

staticprotected

Definition at line 38 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

bool MaxLikelihoodFit::saveNormalizations_ = false

staticprotected

Definition at line 37 of file MaxLikelihoodFit.h.

Referenced by applyOptions(), and runSpecific().

std::string MaxLikelihoodFit::signalPdfNames_ = "shapeSig*"

staticprotected

Definition at line 36 of file MaxLikelihoodFit.h.

Referenced by MaxLikelihoodFit(), and runSpecific().

TTree* MaxLikelihoodFit::t_fit_b_

protected

Definition at line 46 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), and runSpecific().

TTree * MaxLikelihoodFit::t_fit_sb_

protected

Definition at line 46 of file MaxLikelihoodFit.h.

Referenced by createFitResultTrees(), and runSpecific().