ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
ProfileInspector.cxx
Go to the documentation of this file.
1 // @(#)root/roostats:$Id: ProfileInspector.cxx 34109 2010-06-24 15:00:16Z moneta $
2 
3 /*************************************************************************
4  * Project: RooStats *
5  * Package: RooFit/RooStats *
6  * Authors: *
7  * Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke *
8  * Akira Shibata
9  *************************************************************************
10  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
11  * All rights reserved. *
12  * *
13  * For the licensing terms see $ROOTSYS/LICENSE. *
14  * For the list of contributors see $ROOTSYS/README/CREDITS. *
15  *************************************************************************/
16 
17 ////////////////////////////////////////////////////////////////////////////////
18 
19 
21 #include "RooRealVar.h"
22 #include "RooAbsReal.h"
23 #include "RooArgSet.h"
24 #include "RooAbsPdf.h"
25 #include "RooArgSet.h"
26 #include "RooCurve.h"
27 #include "TAxis.h"
28 
29 /// ClassImp for building the THtml documentation of the class
31 
32 using namespace RooStats;
33 using namespace std;
34 
35 ////////////////////////////////////////////////////////////////////////////////
36 
38 {
39 }
40 
41 ////////////////////////////////////////////////////////////////////////////////
42 /// ProfileInspector destructor
43 
45 {
46 }
47 
48 ////////////////////////////////////////////////////////////////////////////////
49 
51 {
52  //
53  // < This tool makes a plot of the conditional maximum likelihood estimate of the nuisance parameter
54  // vs the parameter of interest >
55  //
56  // This enables you to discover if any of the nuisance parameters are behaving strangely
57  // curve is the optional parameters, when used you can specify the points previously scanned
58  // in the process of plotOn or createHistogram.
59  // To do this, you can do the following after the plot has been made:
60 
61  // profile, RooRealVar * poi, RooCurve * curve ){
62  //RooCurve * curve = 0;
63 
64  const RooArgSet* poi_set = config->GetParametersOfInterest();
65  const RooArgSet* nuis_params=config->GetNuisanceParameters();
66  RooAbsPdf* pdf = config->GetPdf();
67 
68 
69  if(!poi_set){
70  cout << "no parameters of interest" << endl;
71  return 0;
72  }
73 
74  if(poi_set->getSize()!=1){
75  cout << "only one parameter of interest is supported currently" << endl;
76  return 0;
77  }
78  RooRealVar* poi = (RooRealVar*) poi_set->first();
79 
80 
81  if(!nuis_params){
82  cout << "no nuisance parameters" << endl;
83  return 0;
84  }
85 
86  if(!pdf){
87  cout << "pdf not set" << endl;
88  return 0;
89  }
90 
91  RooAbsReal* nll = pdf->createNLL(data);
92  RooAbsReal* profile = nll->createProfile(*poi);
93 
94  TList * list = new TList;
95  Int_t curve_N=100;
96  Double_t* curve_x=0;
97 // if(curve){
98 // curve_N=curve->GetN();
99 // curve_x=curve->GetX();
100 // } else {
101  Double_t max = dynamic_cast<RooAbsRealLValue*>(poi)->getMax();
102  Double_t min = dynamic_cast<RooAbsRealLValue*>(poi)->getMin();
103  Double_t step = (max-min)/(curve_N-1);
104  curve_x=new Double_t[curve_N];
105  for(int i=0; i<curve_N; ++i){
106  curve_x[i]=min+step*i;
107  }
108 // }
109 
110  map<string, std::vector<Double_t> > name_val;
111  for(int i=0; i<curve_N; i++){
112  poi->setVal(curve_x[i]);
113  profile->getVal();
114 
115  TIterator* nuis_params_itr=nuis_params->createIterator();
116  TObject* nuis_params_obj;
117  while((nuis_params_obj=nuis_params_itr->Next())){
118  RooRealVar* nuis_param = dynamic_cast<RooRealVar*>(nuis_params_obj);
119  if(nuis_param) {
120  string name = nuis_param->GetName();
121  if(nuis_params->getSize()==0) continue;
122  if(nuis_param && (! nuis_param->isConstant())){
123  if(name_val.find(name)==name_val.end()) name_val[name]=std::vector<Double_t>(curve_N);
124  name_val[name][i]=nuis_param->getVal();
125 
126  if(i==curve_N-1){
127  TGraph* g = new TGraph(curve_N, curve_x, &(name_val[name].front()));
128  g->SetName((name+"_"+string(poi->GetName())+"_profile").c_str());
129  g->GetXaxis()->SetTitle(poi->GetName());
130  g->GetYaxis()->SetTitle(nuis_param->GetName());
131  g->SetTitle("");
132  list->Add(g);
133  }
134  }
135  }
136  }
137  }
138 
139  delete [] curve_x;
140 
141 
142  delete nll;
143  delete profile;
144  return list;
145 }
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Definition: RooAbsPdf.cxx:777
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:52
static Vc_ALWAYS_INLINE int_v min(const int_v &x, const int_v &y)
Definition: vector.h:433
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
Definition: ModelConfig.h:244
int Int_t
Definition: RtypesCore.h:41
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
Definition: ModelConfig.h:250
Iterator abstract base class.
Definition: TIterator.h:32
RooAbsArg * first() const
TIterator * createIterator(Bool_t dir=kIterForward) const
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
A doubly linked list.
Definition: TList.h:47
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:203
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
Bool_t isConstant() const
Definition: RooAbsArg.h:266
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
virtual ~ProfileInspector()
ProfileInspector destructor.
static Vc_ALWAYS_INLINE int_v max(const int_v &x, const int_v &y)
Definition: vector.h:440
#define name(a, b)
Definition: linkTestLib0.cpp:5
Mother of all ROOT objects.
Definition: TObject.h:58
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual void Add(TObject *obj)
Definition: TList.h:81
TList * GetListOfProfilePlots(RooAbsData &data, RooStats::ModelConfig *config)
virtual RooAbsReal * createProfile(const RooArgSet &paramsOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function...
Definition: RooAbsReal.cxx:463
virtual TObject * Next()=0
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
Definition: ModelConfig.h:247
Int_t getSize() const
ClassImp(RooStats::ProfileInspector)
ClassImp for building the THtml documentation of the class.