You are viewing an old version of this page. View the current version.

Compare with Current View Page History

« Previous Version 12 Next »

Abstract

We discuss PSHist - the histogramming package for PSANA project.

Objectives

In PSANA framework we need in package which allows to accumulate data in form of histograms and tuples and save them in file for further analysis. Though it might be based on well-known underlying packages like ROOT, HBOOK, or HippoTuple, we prefer to use an uniform interface with abstract base class, substituting the direct interaction with underlying methods. This intermediate abstract layer provides a flexibility in implementation of new things, for example multithreding in analysis, which algorithm is not yet defined.

Structure and content of the package

Package PSHist contains the base abstract class and methods, which are implemented in derived class RootHist based on CERN's ROOT library. Inheritance between PSHist, RootHist, and Root classes is presented in table

PSHist       RootHist       ROOT
---------------------------------
HManager  <- RootManager <- TFile
H1        <- RootH1      <- TH1
H2        <- RootH2      <- TH2
Profile   <- RootProfile <- TProfile
Tuple     <- RootTuple   <- TTree
Column    <- RootColumn  <- TBranch

Headers

In order to use PSHist one need to include headers:

#include "PSHist/HManager.h"
#include "PSHist/Axis.h"
#include "PSHist/H1.h"
#include "PSHist/H2.h"
#include "PSHist/Profile.h"
#include "PSHist/Tuple.h"
#include "PSHist/Column.h"

ROOT headers which are needed for test purpose only:

#include "RootHist/RootHManager.h" // Is used for initialization
#include "root/TRandom.h"          // For random number generators

Interface

An example of the program interface for this histogram-tuple-management package can be presented as follows. First, the underlying (derived) package needs to be chosen and initialized. For example for ROOT

  PSHist::HManager *hMan = new RootHist::RootHManager("pshist-test.root", "RECREATE"); // for ROOT

In PSANA framework user do not need in this initialization. It will be done automatically using the file name and mode from the list of configuration parameters.

Then later in the program we may create ntuple(s) without knowing what is the underlying histogram-tuple-storage management system,

  PSHist::Tuple *pTuple = hMan->tuple( "TUPLE_N1", "My tuple title^{#alpha}" );

1d histograms

  PSHist::H1 *pH1_1 = hMan->hist1i( "H1_N0001", "My his1d 1 title^{#alpha}", 100, 0., 1.);
  PSHist::H1 *pH1_2 = hMan->hist1f( "H1_N0002", "My his1d 2 title^{#beta}",  axis4 );
  PSHist::H1 *pH1_4 = hMan->hist1d( "H1_N0004", "My his1d 4 title^{#delta}", 10, edges);

2d histograms

  PSHist::H2 *pH2_1 = hMan->hist2i( "H2_N0001", "My his2d 1 title^{#theta}", axis1, axis2 );
  PSHist::H2 *pH2_2 = hMan->hist2f( "H2_N0002", "My his2d 2 title^{#phi}",   axis1, axis3 );

and 1d profile histograms

  // Profile - 1d profile histograms
  PSHist::Profile *pP1_1 = hMan->prof1( "P1_N0001", "My profile 1 title^{#psi}",   60, 0, 1,  0, 1 );
  PSHist::Profile *pP1_2 = hMan->prof1( "P1_N0002", "My profile 2 title^{#iota}",  axis2,     0, 1 );
  PSHist::Profile *pP1_3 = hMan->prof1( "P1_N0003", "My profile 3 title^{#kappa}", 10, edges, 0, 1 );

where we use the class PSHist::Axis for declaration of equal or variable size bins:

  double edges[]={0, 0.05, 0.1, 0.2, 0.25, 0.3, 0.5, 0.55, 0.7, 0.9, 1};
  PSHist::Axis axis1(60,0,1);
  PSHist::Axis axis3(10,edges);

When the ntuple is created, its parameters (Columns) need to be defined once,

  double val;
  float freq;
  typedef struct {float x,y,z;} POINT;
  static POINT point;

  PSHist::Column *pColumn_1 = pTuple->column( &val,  "EBEAM/D" );
                              pTuple->column( &freq, "Freq/F"  );
                              pTuple->column( &point,"x:y:z"   );

      //TuplePar *pBeamEnergy  = nt->parameter("beamEnergy");
      //TuplePar *pBeamCurrent = nt->parameter("beamCurrent");

Then, in data processing stage, for example for each event, we may accumulate data in created histogram and ntuple objects using pointers

            val     = gRandom->Gaus(2.5, 0.1);
            freq    = gRandom->Gaus(1.5, 0.1);
            point.x = gRandom->Gaus(1, 1);
            point.y = gRandom->Gaus(2, 1);
            point.z = gRandom->Gaus(3, 1);

            pTuple -> fill();

pointers

      //pBeamEnergy ->fill(E);
      //pBeamCurrent->fill(I);

or parameter names

      //nt         ->fill("beamEnergy", E);
      //nt         ->fill("beamCurrent", I);
      //nt         ->addRow();

Histograms and profiles can be accumulated using their pointers

            pH1_1 -> fill( gRandom->Rndm(1) ); // Uniform distribution on [0,1]
            pH2_1 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );
            pP1_1 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );

        //pointer -> fill(x,[weight]);
        //pointer -> fill(x,y,[weight]);

where we assume that all parameters like x, y, E, and I and optional weight were earlier defined.

When you are all done, write the data into a file:

  hMan -> write();
  delete hMan;  // close file

In PSANA framework the delete hMan; statement will be done automatically.

Code examples

Accumulation of data

//-----------------
// C/C++ Headers --
//-----------------

#include <string>
#include <iostream>

using std::cout;
using std::endl;
//----------------------
// Headers
//----------------------

#include "RootHist/RootHManager.h"
#include "PSHist/HManager.h"
#include "PSHist/Axis.h"
#include "PSHist/H1.h"
#include "PSHist/H2.h"
#include "PSHist/Profile.h"
#include "PSHist/Tuple.h"
#include "PSHist/Column.h"

#include "root/TRandom.h"

int main ()
{
  cout << "Start main()" << endl;

  // Non-equal bin edges:
  double edges[]={0, 0.05, 0.1, 0.2, 0.25, 0.3, 0.5, 0.55, 0.7, 0.9, 1};

  // Axis
  PSHist::Axis axis1(60,0,1);
  PSHist::Axis axis2(40,0,1);
  PSHist::Axis axis3(10,edges);
  PSHist::Axis axis4(100,0,1);

  PSHist::HManager *hMan = new RootHist::RootHManager("pshist-test.root", "RECREATE");

  // H1 - 1d histograms
  PSHist::H1 *pH1_1 = hMan->hist1i( "H1_N0001", "My his1d 1 title^{#alpha}", 100, 0., 1.);
  PSHist::H1 *pH1_2 = hMan->hist1f( "H1_N0002", "My his1d 2 title^{#beta}",  axis4 );
  PSHist::H1 *pH1_3 = hMan->hist1d( "H1_N0003", "My his1d 3 title^{#gamma}", 100, 0., 1.01);
  PSHist::H1 *pH1_4 = hMan->hist1d( "H1_N0004", "My his1d 4 title^{#delta}", 10, edges);

  // H2 - 2d histograms
  PSHist::H2 *pH2_1 = hMan->hist2i( "H2_N0001", "My his2d 1 title^{#theta}",  axis1, axis2 );
  PSHist::H2 *pH2_2 = hMan->hist2f( "H2_N0002", "My his2d 2 title^{#phi}",    axis1, axis3 );
  PSHist::H2 *pH2_3 = hMan->hist2d( "H2_N0003", "My his2d 3 title^{#eta}",    axis2, axis2 );
  PSHist::H2 *pH2_4 = hMan->hist2d( "H2_N0004", "My his2d 4 title^{#lambda}", axis3, axis2 );

  // Profile - 1d profile histograms
  PSHist::Profile *pP1_1 = hMan->prof1( "P1_N0001", "My profile 1 title^{#psi}",     60, 0, 1,  0, 1 );
  PSHist::Profile *pP1_2 = hMan->prof1( "P1_N0002", "My profile 2 title^{#iota}",    axis2,     0, 1 );
  PSHist::Profile *pP1_3 = hMan->prof1( "P1_N0003", "My profile 3 title^{#kappa}",   10, edges, 0, 1 );
  PSHist::Profile *pP1_4 = hMan->prof1( "P1_N0004", "My profile 4 title^{#lambda}",  axis3,     0, 1 );

  // Tuple
  PSHist::Tuple *pTuple = hMan->tuple( "TUPLE_N1", "My tuple title^{#alpha}" );

  // Column (parameter(s)) for the tuple; use the ROOT-style constructor
  double val;
  float freq;
  typedef struct {float x,y,z;} POINT;
  static POINT point;
  /*
  PSHist::Column *pColumn_1 = pTuple->column( "C_N0001", &val,  "EBEAM/D" );
                              pTuple->column( "C_N0002", &freq, "Freq/F"  );
                              pTuple->column( "point",   &point,"x:y:z"   );
  */
  PSHist::Column *pColumn_1 = pTuple->column( &val,  "EBEAM/D" );
                              pTuple->column( &freq, "Freq/F"  );
                              pTuple->column( &point,"x:y:z"   );

                  pColumn_1->print(cout);

  cout << "Fill histograms" << endl;
        for (int i=0; i<10000; i++)
          {
            pH1_1 -> fill( gRandom->Rndm(1)        ); // Uniform distribution on [0,1]
            pH1_2 -> fill( gRandom->Gaus(0.5, 0.1) ); // Gaussian for mu and sigma
            pH1_3 -> fill( double(0.0001*i), double(i) );
            pH1_4 -> fill( gRandom->Rndm(1)        );

            pH2_1 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );
            pH2_2 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );
            pH2_3 -> fill( gRandom->Gaus(0.4, 0.1), gRandom->Gaus(0.6, 0.2) );
            pH2_4 -> fill( gRandom->Gaus(0.3, 0.2), gRandom->Gaus(0.7, 0.1) );

            pP1_1 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );
            pP1_2 -> fill( gRandom->Gaus(0.5, 0.1), gRandom->Rndm(1) );
            pP1_3 -> fill( gRandom->Gaus(0.4, 0.1), gRandom->Gaus(0.6, 0.2) );
            pP1_4 -> fill( gRandom->Gaus(0.3, 0.2), gRandom->Gaus(0.7, 0.1) );

            val  = gRandom->Gaus(2.5, 0.1);
            freq = gRandom->Gaus(1.5, 0.1);

            point.x = gRandom->Gaus(1, 1);
            point.y = gRandom->Gaus(2, 1);
            point.z = gRandom->Gaus(3, 1);

            pTuple -> fill();
          }

  hMan -> write();
  delete hMan;  // close file
}

Browse the ROOT file

In order to browser data stored in file we use a ROOT script proc-pshist.C:

//void proc(int Nplot=1)
{
  // Settings for good style
   gStyle -> SetPadColor(3);
   gStyle -> SetPadBorderSize(0);
   gStyle -> SetPadBorderMode(0);
   gStyle -> SetTitleXSize(0.05); // set size of axes titles
   gStyle -> SetTitleYSize(0.05);
   gStyle -> SetTitleH(0.1); // set size of the title in top box

  TFile *f = new TFile("pshist-test.root");
         f -> ls();

  //------------------------------------

  c1 = new TCanvas("c1","",0,0,800,800);
  c1->Divide(2,2);

  c1->cd(1); H1_N0001->Draw();
  c1->cd(2); H1_N0002->Draw();
  c1->cd(3); H1_N0003->Draw();
  c1->cd(4); H1_N0004->Draw();

  gPad -> Update();

  //------------------------------------

  c2 = new TCanvas("c2","",200,0,800,800);
  c2->Divide(2,2);

  c2->cd(1); H2_N0001->Draw("box");
  c2->cd(2); H2_N0002->Draw("COLORZ");
  c2->cd(3); H2_N0003->Draw("box");
  c2->cd(4); H2_N0004->Draw("COLORZ");

  gPad -> Update();

  //------------------------------------

  c3 = new TCanvas("c3","",400,0,800,800);
  c3->Divide(2,2);

  c3->cd(1); P1_N0001->Draw();
  c3->cd(2); P1_N0002->Draw();
  c3->cd(3); P1_N0003->Draw();
  c3->cd(4); P1_N0004->Draw();

  gPad -> Update();

  //------------------------------------

  c4 = new TCanvas("c4","",600,0,800,800);
  c4->Divide(2,2);

             TUPLE_N1->Print();
  c4->cd(1); TUPLE_N1->Draw("EBEAM");
  c4->cd(2); TUPLE_N1->Draw("Freq");
  c4->cd(3); TUPLE_N1->Draw("x");
  c4->cd(4); TUPLE_N1->Draw("y");

  gPad -> Update();

  //------------------------------------

  c1 -> Print("plot1.gif"); // works for gif, eps, etc.
  c2 -> Print("plot2.gif");
  c3 -> Print("plot3.gif");
  c4 -> Print("plot4.gif");

  //------------------------------------

  cout << "Sleep for 10 sec..." << endl;
  gSystem->Sleep(10*1000);
  cout << "Wake up, and exit root." << endl;

  f -> Close();
}

To run this script use command

root -q -f proc-pshist.C

This script draws four canvas and saves them in files plot1.gif-plot4.gif,

  • No labels