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

Compare with Current View Page History

« Previous Version 2 Next »

Abstract

We discuss a PSHist - a 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.

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/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"

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 choisen and instantiated. For example for ROOT

      HManager *hMan = new RootHist("my-data-file.root"); // for ROOT

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

      Tuple *nt = hMan->ntuple("EXP Data");

and histograms

     H1 *pHis1i = hMan->his1i("His1 int    title",100,0.,1.);
     H1 *pHis1f = hMan->his1f("His1 float  title",100,0.,1.);
     H1 *pHis1d = hMan->his1d("His1 double title",100,0.,1.);
     H2 *pHis2d = hMan->his2d("His2 double title",100,0.,1.,100,0.,1.);

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

      TuplePar *pBeamEnergy  = nt->parameter("beamEnergy");  // minValue, maxValue, dtype ?
      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

      pBeamEnergy ->fill(E);
      pBeamCurrent->fill(I);
      nt          ->addRow();

or parameter names

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

Histograms can be accumulated using their pointers

      pHis1i     ->fill(x,[weight]);
      ...
      pHis2      ->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;

Structure and content of the package

Package PSHist contains the base abstract class and methods, which have to be re-implemented in derived classes for HBOOK, ROOT, HippoTuple, etc., i.e. in packages like HBookHist, RootHist, HippoHist, etc.

Inheritance

Mapping between PSHist and RootHist and Root classes:

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

Code examples

Accumulation of data

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

#include <string>
#include <iostream>
//#include <sstream>

//#include <stdio.h> // for C style printf

//----------------------
// Base Class Headers --
//----------------------

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

//-------------------------------
// Collaborating Class Headers --
//-------------------------------

#include "root/TRandom.h"

//-----------------------------------------------------------------------
// Local Macros, Typedefs, Structures, Unions and Forward Declarations --
//-----------------------------------------------------------------------

using std::cout;
using std::endl;

//using namespace PSTime;

int main ()
{
  cout << "Stars ex_PSHist_1 : main()" << endl;

  double edges[]={0, 0.05, 0.1, 0.2, 0.25, 0.3, 0.5, 0.55, 0.7, 0.9, 1};

  cout << "sizeof(edges) =" << sizeof(edges) << endl;


  // 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

// To run this script use command: root -q -f proc-pshist.C


//void proc(int Nplot=1)
{
  // Settings for pads

   gStyle->SetPadColor(3);
   gStyle->SetPadBorderSize(0);
   gStyle->SetPadBorderMode(0);

   gStyle -> SetTitleXSize(0.05); // set size of axes titles
   gStyle -> SetTitleYSize(0.05);

   //gStyle->SetTitleW(0.35);
   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();

}

  • No labels