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(); }