%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooDataHist.h"
#include "RooGaussian.h"
#include "RooCategory.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "TFile.h"
using namespace RooFit;

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

RooDataSet d("d", "d", RooArgSet(x, y, c));

Int_t i;
for (i = 0; i < 1000; i++) {
   x = i / 50 - 10;
   y = sqrt(1.0 * i);
   c.setLabel((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;

const RooArgSet *row = d.get();
row->Print("v");
cout << endl;

d.get(900)->Print("v");
cout << endl;

cout << endl << ">> d1 has only columns x,c" << endl;
std::unique_ptr<RooAbsData> d1{d.reduce({x, c})};
d1->Print("v");

cout << endl << ">> d2 has only column y" << endl;
std::unique_ptr<RooAbsData> d2{d.reduce({y})};
d2->Print("v");

cout << endl << ">> d3 has only the points with y>5.17" << endl;
std::unique_ptr<RooAbsData> 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;
std::unique_ptr<RooAbsData> d4{d.reduce({x, c}, "y>5.17")};
d4->Print("v");

cout << endl << ">> merge d2(y) with d1(x,c) to form d1(x,c,y)" << endl;
static_cast<RooDataSet&>(*d1).merge(&static_cast<RooDataSet&>(*d2));
d1->Print("v");

cout << endl << ">> append data points of d3 to d1" << endl;
static_cast<RooDataSet&>(*d1).append(static_cast<RooDataSet&>(*d3));
d1->Print("v");

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;

x.setBins(10);
y.setBins(10);
RooDataHist dh("dh", "binned version of d", RooArgSet(x, y), d);
dh.Print("v");

RooPlot *yframe = y.frame(Bins(10), Title("Operations on binned datasets"));
dh.plotOn(yframe); // plot projection of 2D binned data on y

cout << ">> number of bins in dh   : " << dh.numEntries() << endl;
cout << ">> sum of weights in dh   : " << dh.sum(false) << endl;
cout << ">> integral over histogram: " << dh.sum(true) << endl; // accounts for bin volume

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

cout << ">> Creating 1-dimensional projection on y of dh for bins with x>0" << endl;
std::unique_ptr<RooAbsData> dh2{dh.reduce(y, "x>0")};
dh2->Print("v");

dh2->plotOn(yframe, LineColor(kRed), MarkerColor(kRed));

cout << endl << ">> Persisting d via ROOT I/O" << endl;
TFile f("rf402_datahandling.root", "RECREATE");
d.Write();
f.ls();

new TCanvas("rf402_datahandling", "rf402_datahandling", 600, 600);
gPad->SetLeftMargin(0.15);
yframe->GetYaxis()->SetTitleOffset(1.4);
yframe->Draw();

%jsroot on
gROOT->GetListOfCanvases()->Draw()