void convolution_hist () {
  // macro to create a MC dataset from an exponential function convoluted 
  // with a smeared distribution for the tau constant
  // steps
  // 1) load the histogram containg the possible tau values
  // 2) build a RooHistPdf from a RooHistDataSet based on the histogram
  // 3) create an exponential pdf
  // 4) convolve the exp pdf and the smearing function (use RooFFTConvPdf)
  // 5) create a MC sample based on the convolved Pdf
  //
  gSystem->Load("libRooFit.so");
  using namespace RooFit;

  TCanvas *c1 = new TCanvas("c1","",1024,768);
  c1->cd();

  // 1) load histogram
  TFile* f= new TFile("tests_140224-fei4-250um-fl3.8e15-extract_traps-V400-taue-in-ns.root","READ");
  f->cd();


  TH1D* hTauE = (TH1D*)f->Get("hTauE");
  if ( !hTauE) {
    std::cout << "hTauE not found.\Exiting...\n";
    exit(1);
  }

  // determining tau range
  int nbins = hTauE->GetNbinsX();

  double xmin = hTauE->GetBinLowEdge(1);
  double xmax = hTauE->GetBinLowEdge(nbins+1);
  double mean = hTauE->GetMean();

  std::cout << "mean = " << mean << "\n";
  std::cout << "xmin = " << xmin << "\n";
  std::cout << "xmax = " << xmax << "\n";

  // the tau parameter
  RooRealVar tau("tau","#tau [ns]",mean,xmin,xmax);
  // the data set
  RooDataHist* tauEDH = new RooDataHist("tauEDH","tauEDH",tau,Import(*hTauE));
 
  // 2) build the RooHistPdf
  RooHistPdf tauHPDF("tauHPDF","tauHPDF",tau,*tauEDH,0);
  // check the result
  RooPlot* frame1 = tau.frame(Title("Tau histogram pdf"));
  tauEDH->plotOn(frame1);
  tauHPDF.plotOn(frame1);
  frame1->Draw();

  // 3) The exponential pdf
  RooRealVar t("t","t [ns]",0,fabs(xmax)*5);
  t.Print();
  tau.Print();
  RooExponential expPDF("expPDF","expPDF",t,tau);
  RooPlot* frame2 = t.frame(Title("Trapping time"));
  expPDF.plotOn(frame2);  

  TCanvas *c2 = new TCanvas("c2","",1024,768);
  c2->cd();
  frame2->Draw();

  // 4) convolve the exp pdf and the smearing function (use RooFFTConvPdf)
  // inspired by rf211_paramconv.C
  // Convolution in tau parameter trapConv = expPDF(t,tau) (x) tauHPDF(tau)
  t.setBins(1000,"cache");
  tau.setBins(50,"cache");
  RooFFTConvPdf trapConv("trapConv","trapConv",tau,expPDF,tauHPDF);

  //trapConv.Print();

// Configure convolution to construct a 2-D cache in (t,tau)
// rather than a 1-d cache in tau that needs to be recalculated
// for each value of t
  trapConv.setCacheObservables(t) ;
  trapConv.setBufferFraction(1.0) ;

  // Integrate trapConv over tau; projModel = Int trapConv dtau
  RooAbsPdf* projtrapConv = trapConv.createProjection(tau) ;

  // Generate 1000 toy events
  RooDataHist* d = projtrapConv->generateBinned(t,1000) ;

  // Fit p.d.f. to toy data
  projtrapConv->fitTo(*d,Verbose()) ;

  // Plot data and fitted p.d.f.
  RooPlot* frame3 = t.frame(Bins(25)) ;
  d->plotOn(frame3) ;
  projtrapConv->plotOn(frame3) ;

  TCanvas *c3 = new TCanvas("c3","",1024,768);
  c3->cd();
  frame3->Draw();

  // Make 2d histogram of trapConv(x;mean)
  TH1* hh = trapConv.createHistogram("hh",t,Binning(50),YVar(tau,Binning(50)),ConditionalObservables(tau)) ;
  hh->SetTitle("histogram of trapConv(t|tau)") ;
  hh->SetLineColor(kBlue) ;

  // Draw frame on canvas
  TCanvas* c = new TCanvas("rf211_paramconv","rf211_paramconv",800,400) ;
  c->Divide(2) ;
  c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame3->GetYaxis()->SetTitleOffset(1.4) ; frame3->Draw() ;
  c->cd(2) ; gPad->SetLeftMargin(0.20) ; hh->GetZaxis()->SetTitleOffset(2.5) ; hh->Draw("surf") ;

  // closing the file 
//  f->Close();
}