ROOT logo

From $ROOTSYS/tutorials/roofit/rf705_linearmorph.C

//////////////////////////////////////////////////////////////////////////
//
// 'SPECIAL PDFS' RooFit tutorial macro #705
// 
// Linear interpolation between p.d.f shapes using the 'Alex Read' algorithm
//
//
//
// 07/2008 - Wouter Verkerke 
// 
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooPolynomial.h"
#include "RooIntegralMorph.h"
#include "RooNLLVar.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "TH1.h"
using namespace RooFit ;


void rf705_linearmorph()
{
  // C r e a t e   e n d   p o i n t   p d f   s h a p e s
  // ------------------------------------------------------

  // Observable
  RooRealVar x("x","x",-20,20) ;

  // Lower end point shape: a Gaussian
  RooRealVar g1mean("g1mean","g1mean",-10) ;
  RooGaussian g1("g1","g1",x,g1mean,RooConst(2)) ;

  // Upper end point shape: a Polynomial
  RooPolynomial g2("g2","g2",x,RooArgSet(RooConst(-0.03),RooConst(-0.001))) ;


 
  // C r e a t e   i n t e r p o l a t i n g   p d f 
  // -----------------------------------------------

  // Create interpolation variable
  RooRealVar alpha("alpha","alpha",0,1.0) ;

  // Specify sampling density on observable and interpolation variable
  x.setBins(1000,"cache") ;
  alpha.setBins(50,"cache") ;

  // Construct interpolating pdf in (x,a) represent g1(x) at a=a_min
  // and g2(x) at a=a_max
  RooIntegralMorph lmorph("lmorph","lmorph",g1,g2,x,alpha) ;



  // P l o t   i n t e r p o l a t i n g   p d f   a t   v a r i o u s   a l p h a 
  // -----------------------------------------------------------------------------

  // Show end points as blue curves
  RooPlot* frame1 = x.frame() ;
  g1.plotOn(frame1) ;
  g2.plotOn(frame1) ;

  // Show interpolated shapes in red
  alpha.setVal(0.125) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.25) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.375) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.50) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.625) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.75) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.875) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;
  alpha.setVal(0.95) ;
  lmorph.plotOn(frame1,LineColor(kRed)) ;



  // S h o w   2 D   d i s t r i b u t i o n   o f   p d f ( x , a l p h a ) 
  // -----------------------------------------------------------------------
  
  // Create 2D histogram
  TH1* hh = lmorph.createHistogram("hh",x,Binning(40),YVar(alpha,Binning(40))) ;
  hh->SetLineColor(kBlue) ;


  // F i t   p d f   t o   d a t a s e t   w i t h   a l p h a = 0 . 8 
  // -----------------------------------------------------------------

  // Generate a toy dataset at alpha = 0.8
  alpha=0.8 ;
  RooDataSet* data = lmorph.generate(x,1000) ;

  // Fit pdf to toy data
  lmorph.setCacheAlpha(kTRUE) ;
  lmorph.fitTo(*data,Verbose(kTRUE)) ;

  // Plot fitted pdf and data overlaid
  RooPlot* frame2 = x.frame(Bins(100)) ;
  data->plotOn(frame2) ;
  lmorph.plotOn(frame2) ; 


  // S c a n   - l o g ( L )   v s   a l p h a
  // -----------------------------------------

  // Show scan -log(L) of dataset w.r.t alpha
  RooPlot* frame3 = alpha.frame(Bins(100),Range(0.1,0.9)) ;
  
  // Make 2D pdf of histogram  
  RooNLLVar nll("nll","nll",lmorph,*data) ;  
  nll.plotOn(frame3,ShiftToZero()) ;    

  lmorph.setCacheAlpha(kFALSE) ;



  TCanvas* c = new TCanvas("rf705_linearmorph","rf705_linearmorph",800,800) ;
  c->Divide(2,2) ;
  c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame1->GetYaxis()->SetTitleOffset(1.6) ; frame1->Draw() ;
  c->cd(2) ; gPad->SetLeftMargin(0.20) ; hh->GetZaxis()->SetTitleOffset(2.5) ; hh->Draw("surf") ;
  c->cd(3) ; gPad->SetLeftMargin(0.15) ; frame3->GetYaxis()->SetTitleOffset(1.4) ; frame3->Draw() ;
  c->cd(4) ; gPad->SetLeftMargin(0.15) ; frame2->GetYaxis()->SetTitleOffset(1.4) ; frame2->Draw() ;
  
  
  return ;

}
 rf705_linearmorph.C:1
 rf705_linearmorph.C:2
 rf705_linearmorph.C:3
 rf705_linearmorph.C:4
 rf705_linearmorph.C:5
 rf705_linearmorph.C:6
 rf705_linearmorph.C:7
 rf705_linearmorph.C:8
 rf705_linearmorph.C:9
 rf705_linearmorph.C:10
 rf705_linearmorph.C:11
 rf705_linearmorph.C:12
 rf705_linearmorph.C:13
 rf705_linearmorph.C:14
 rf705_linearmorph.C:15
 rf705_linearmorph.C:16
 rf705_linearmorph.C:17
 rf705_linearmorph.C:18
 rf705_linearmorph.C:19
 rf705_linearmorph.C:20
 rf705_linearmorph.C:21
 rf705_linearmorph.C:22
 rf705_linearmorph.C:23
 rf705_linearmorph.C:24
 rf705_linearmorph.C:25
 rf705_linearmorph.C:26
 rf705_linearmorph.C:27
 rf705_linearmorph.C:28
 rf705_linearmorph.C:29
 rf705_linearmorph.C:30
 rf705_linearmorph.C:31
 rf705_linearmorph.C:32
 rf705_linearmorph.C:33
 rf705_linearmorph.C:34
 rf705_linearmorph.C:35
 rf705_linearmorph.C:36
 rf705_linearmorph.C:37
 rf705_linearmorph.C:38
 rf705_linearmorph.C:39
 rf705_linearmorph.C:40
 rf705_linearmorph.C:41
 rf705_linearmorph.C:42
 rf705_linearmorph.C:43
 rf705_linearmorph.C:44
 rf705_linearmorph.C:45
 rf705_linearmorph.C:46
 rf705_linearmorph.C:47
 rf705_linearmorph.C:48
 rf705_linearmorph.C:49
 rf705_linearmorph.C:50
 rf705_linearmorph.C:51
 rf705_linearmorph.C:52
 rf705_linearmorph.C:53
 rf705_linearmorph.C:54
 rf705_linearmorph.C:55
 rf705_linearmorph.C:56
 rf705_linearmorph.C:57
 rf705_linearmorph.C:58
 rf705_linearmorph.C:59
 rf705_linearmorph.C:60
 rf705_linearmorph.C:61
 rf705_linearmorph.C:62
 rf705_linearmorph.C:63
 rf705_linearmorph.C:64
 rf705_linearmorph.C:65
 rf705_linearmorph.C:66
 rf705_linearmorph.C:67
 rf705_linearmorph.C:68
 rf705_linearmorph.C:69
 rf705_linearmorph.C:70
 rf705_linearmorph.C:71
 rf705_linearmorph.C:72
 rf705_linearmorph.C:73
 rf705_linearmorph.C:74
 rf705_linearmorph.C:75
 rf705_linearmorph.C:76
 rf705_linearmorph.C:77
 rf705_linearmorph.C:78
 rf705_linearmorph.C:79
 rf705_linearmorph.C:80
 rf705_linearmorph.C:81
 rf705_linearmorph.C:82
 rf705_linearmorph.C:83
 rf705_linearmorph.C:84
 rf705_linearmorph.C:85
 rf705_linearmorph.C:86
 rf705_linearmorph.C:87
 rf705_linearmorph.C:88
 rf705_linearmorph.C:89
 rf705_linearmorph.C:90
 rf705_linearmorph.C:91
 rf705_linearmorph.C:92
 rf705_linearmorph.C:93
 rf705_linearmorph.C:94
 rf705_linearmorph.C:95
 rf705_linearmorph.C:96
 rf705_linearmorph.C:97
 rf705_linearmorph.C:98
 rf705_linearmorph.C:99
 rf705_linearmorph.C:100
 rf705_linearmorph.C:101
 rf705_linearmorph.C:102
 rf705_linearmorph.C:103
 rf705_linearmorph.C:104
 rf705_linearmorph.C:105
 rf705_linearmorph.C:106
 rf705_linearmorph.C:107
 rf705_linearmorph.C:108
 rf705_linearmorph.C:109
 rf705_linearmorph.C:110
 rf705_linearmorph.C:111
 rf705_linearmorph.C:112
 rf705_linearmorph.C:113
 rf705_linearmorph.C:114
 rf705_linearmorph.C:115
 rf705_linearmorph.C:116
 rf705_linearmorph.C:117
 rf705_linearmorph.C:118
 rf705_linearmorph.C:119
 rf705_linearmorph.C:120
 rf705_linearmorph.C:121
 rf705_linearmorph.C:122
 rf705_linearmorph.C:123
 rf705_linearmorph.C:124
 rf705_linearmorph.C:125
 rf705_linearmorph.C:126
 rf705_linearmorph.C:127
 rf705_linearmorph.C:128
 rf705_linearmorph.C:129
 rf705_linearmorph.C:130
 rf705_linearmorph.C:131
 rf705_linearmorph.C:132
 rf705_linearmorph.C:133
 rf705_linearmorph.C:134
 rf705_linearmorph.C:135
 rf705_linearmorph.C:136
 rf705_linearmorph.C:137
 rf705_linearmorph.C:138
 rf705_linearmorph.C:139
 rf705_linearmorph.C:140
 rf705_linearmorph.C:141