JTools/JHistogram2D.cc File Reference

#include <string>
#include <iostream>
#include <iomanip>
#include "TROOT.h"
#include "TFile.h"
#include "TH2D.h"
#include "TMath.h"
#include "TRandom.h"
#include "JTools/JHistogram1D_t.hh"
#include "JTools/JHistogramMap_t.hh"
#include "JTools/JCollection.hh"
#include "JTools/JMultiGrid.hh"
#include "JTools/JMultiSet.hh"
#include "JTools/JQuadrature.hh"
#include "JTools/JFunction1D_t.hh"
#include "JTools/JFunctionalMap_t.hh"
#include "JTools/JToolsToolkit.hh"
#include "JTools/JMultiPDF.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions
int	main (int argc, char **argv)

Detailed Description

Example program to test JTOOLS::JHistogram2D.

Author

mdejong

Definition in file JTools/JHistogram2D.cc.

Function Documentation

main()

int main	(	int	argc,
		char **	argv
	)

Definition at line 84 of file JTools/JHistogram2D.cc.

{
  using namespace std;
 
  string  outputFile;
  int     numberOfEvents;
  int     numberOfBins;
  bool    quadrature;
  bool    subtract;
  int     debug;
 
  try {
 
    JParser<> zap("Example program to test 2D histogram.");
 
    zap['o'] = make_field(outputFile)     = "histogram.root";
    zap['n'] = make_field(numberOfEvents);
    zap['N'] = make_field(numberOfBins);
    zap['Q'] = make_field(quadrature);
    zap['D'] = make_field(subtract);
    zap['d'] = make_field(debug)          = 3;
 
    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }
 
 
  if (numberOfEvents <= 0) {
    FATAL("No events." << endl);
  }
 
  using namespace JPP;
 
  const Int_t    nx    = 1000;
  const Double_t xmin  = -3.0;
  const Double_t xmax  = +3.0;
 
  typedef JMultiHistogram<JHistogram1D_t,
                          JMapList<JHistogramMap_t> >  JHistogram2D_t;
 
  JHistogram2D_t histogram;
  
 
  // abscissa
 
  if (!quadrature)
    configure(histogram, make_grid(numberOfBins + 1, xmin, xmax));
  else
    configure(histogram, make_set(JGaussHermite(numberOfBins + 1)));
 
  
  // fill
 
  for (int i = 0; i != numberOfEvents; ++i) {
 
    const double x = gRandom->Gaus(0.0, 1.0);
    const double y = gRandom->Gaus(0.0, 1.0);
 
    histogram.fill(x, y, 1.0);
  }
  
  
  histogram.div((double) numberOfEvents);
 
  
  // interpolation based on function values
 
  JMultiPDF<JPolint2Function1D_t,
            JMapList<JPolint2FunctionalMap> >         f2(histogram);
 
 
  // interpolation based on integral values
 
  accumulate(histogram);
  
  JMultiFunction<JPolint3Function1H_t,
                 JMapList<JPolint3FunctionalMapH> >   F2;
 
  copy(histogram, F2);
 
  F2.compile();
 
 
  TFile out(outputFile.c_str(), "recreate");
 
  TH2D h0("h0", NULL, nx, xmin, xmax, nx, xmin, xmax);
  TH2D h1("h1", NULL, nx, xmin, xmax, nx, xmin, xmax);
  TH2D h2("h2", NULL, nx, xmin, xmax, nx, xmin, xmax);
 
  
  for (int i = 1; i <= h0.GetXaxis()->GetNbins(); ++i) {
    for (int j = 1; j <= h0.GetYaxis()->GetNbins(); ++j) {
      
      const Double_t x = h0.GetXaxis()->GetBinCenter(i);
      const Double_t y = h0.GetYaxis()->GetBinCenter(j);
      
      h0.SetBinContent(i, j, g2(x,y));
      
      try {
        h1.SetBinContent(i, j, F2(x,y).fp.fp);
        h2.SetBinContent(i, j, f2(x,y));
      }
      catch(const exception& error) {
        //ERROR(error.what() << endl);
      }
    }
  }
 
  if (subtract) {
    h1.Add(&h0, -1.0);
    h2.Add(&h0, -1.0);
  }
 
  out.Write();
  out.Close();
}