Docs | All | Real | Category | PDF | DataSet | Plot | Container | Misc | Aux | User

RooFit Toolkit for Data Modeling

Dataset operations

// Dataset operations
intro8()
{
  // Binned (RooDataHist) and unbinned datasets (RooDataSet) share
  // many properties and inherit from a common abstract base class
  // (RooAbsData), that provides an interface for all operations
  // that can be performed regardless of the data format

  RooRealVar  x("x","x",-10,10) ;
  RooRealVar  y("y","y", 0, 40) ;
  RooCategory c("c","c") ;
  c.defineType("Plus",+1) ;
  c.defineType("Minus",-1) ;

  // *** Unbinned datasets ***

  // RooDataSet is an unbinned dataset (a collection of points in N-dimensional space)
  RooDataSet d("d","d",RooArgSet(x,y,c)) ;

  // Unlike RooAbsArgs (RooAbsPdf,RooFormulaVar,....) datasets are not attached to 
  // the variables they are constructed from. Instead they are attached to an internal 
  // clone of the supplied set of arguments

  // Fill d with dummy values
  Int_t i ;
  for (i=0 ; i<1000 ; i++) {
    x = i/50 - 10 ;
    y = sqrt(1.0*i) ;
    c = (i%2)?"Plus":"Minus" ;

    // We must explicitly refer to x,y,c here to pass the values because
    // d is not linked to them (as explained above)
    d.add(RooArgSet(x,y,c)) ;
  }
  d.Print("v") ;
  cout << endl ;

  // The get() function returns a pointer to the internal copy of the RooArgSet(x,y,c)
  // supplied in the constructor
  RooArgSet* row = d.get() ;
  row->Print("v") ;
  cout << endl ;

  // Get with an argument loads a specific data point in row and returns
  // a pointer to row argset. get() always returns the same pointer, unless
  // an invalid row number is specified. In that case a null ptr is returned
  d.get(900)->Print("v") ;
  cout << endl ;

  // *** Reducing / Appending / Merging ***

  // The reduce() function returns a new dataset which is a subset of the original
  cout << endl << ">> d1 has only columns x,c" << endl ;
  RooDataSet* d1 = d.reduce(RooArgSet(x,c)) ; 
  d1->Print("v") ;

  cout << endl << ">> d2 has only column y" << endl ;
  RooDataSet* d2 = d.reduce(RooArgSet(y)) ;   
  d2->Print("v") ;

  cout << endl << ">> d3 has only the points with y>5.17" << endl ;
  RooDataSet* d3 = d.reduce("y>5.17") ; 
  d3->Print("v") ;

  cout << endl << ">> d4 has only columns x,c for data points with y>5.17" << endl ;
  RooDataSet* d4 = d.reduce(RooArgSet(x,c),"y>5.17") ; 
  d4->Print("v") ;

  // The merge() function adds two data set column-wise
  cout << endl << ">> merge d2(y) with d1(x,c) to form d1(x,c,y)" << endl ;
  d1->merge(d2) ; 
  d1->Print("v") ;

  // The append() function addes two datasets row-wise
  cout << endl << ">> append data points of d3 to d1" << endl ;
  d1->append(*d3) ;
  d1->Print("v") ;  

  // *** Binned datasets ***

  // A binned dataset can be constructed empty, from an unbinned dataset, or
  // from a ROOT native histogram (TH1,2,3)

  cout << ">> construct dh (binned) from d(unbinned) but only take the x and y dimensions," << endl 
       << ">> the category 'c' will be projected in the filling process" << endl ;
    
  // The binning of real variables (like x,y) is done using their fit range
  //'get/setFitRange()' and number of specified fit bins 'get/setFitBins()'.
  // Category dimensions of binned datasets get one bin per defined category state
  x.setFitBins(10) ;
  y.setFitBins(10) ;
  RooDataHist dh("dh","binned version of d",RooArgSet(x,y),d) ;
  dh.Print("v") ;

  RooPlot* yframe = y.frame(10) ;
  dh.plotOn(yframe) ; // plot projection of 2D binned data on y
  yframe->Draw() ;

  // Examine the statistics of a binned dataset
  cout << ">> number of bins in dh   : " << dh.numEntries() << endl ;
  cout << ">> sum of weights in dh   : " << dh.sum(kFALSE)  << endl ;
  cout << ">> integral over histogram: " << dh.sum(kTRUE)   << endl ; // accounts for bin volume

  // Locate a bin from a set of coordinates and retrieve its properties
  x = 0.3 ;  y = 20.5 ;
  cout << ">> retrieving the properties of the bin enclosing coordinate (x,y) = (0.3,20.5) " << endl ;
  cout << " bin center:" << endl ;
  dh.get(RooArgSet(x,y))->Print("v") ; // load bin center coordinates in internal buffer
  cout << " weight = " << dh.weight() << endl ; // return weight of last loaded coordinates

  // Reduce the 2-dimensional binned dataset to a 1-dimensional binned dataset
  //
  // All reduce() methods are interfaced in RooAbsData. All reduction techniques
  // demonstrated on unbinned datasets can be applied to binned datasets as well.
  cout << ">> Creating 1-dimensional projection on y of dh for bins with x>0" << endl ;
  RooDataHist* dh2 = dh.reduce(y,"x>0") ;
  dh2->Print("v") ;

  // Add dh2 to yframe and redraw
  dh2->plotOn(yframe,LineColor(kRed),MarkerColor(kRed)) ;
  yframe->Draw() ;

  // *** Saving and loading from file ***
  
  // Datasets can be persisted with ROOT I/O
  cout << endl << ">> Persisting d via ROOT I/O" << endl ;
  TFile f("intro8.root","RECREATE") ;
  d.Write() ;
  f.ls() ;

  // To read back in future session:
  // > TFile f("intro8.root") ;
  // > RooDataSet* d = (RooDataSet*) f.FindObject("d") ;
  
}
Last CVS Update: Top
Copyright © 2000-2005 University of California, Stanford University

Page maintained by Wouter Verkerke and David Kirkby

SourceForge.net Logo