ROOT
master
Reference Guide
Loading...
Searching...
No Matches
hsimple.py
Go to the documentation of this file.
1
## \file
2
## \ingroup Tutorials
3
## \notebook -js
4
## This program creates :
5
## - a one dimensional histogram
6
## - a two dimensional histogram
7
## - a profile histogram
8
## - a memory-resident ntuple
9
##
10
## These objects are filled with some random numbers and saved on a file.
11
##
12
## \macro_image
13
## \macro_code
14
##
15
## \author Wim Lavrijsen, Enric Tejedor
16
## \author Vincenzo Eduardo Padulano (CERN), 09.2025
17
18
import
numpy
19
from
ROOT
import
TH1F, TH2F, TCanvas, TFile, TNtuple, TProfile, gBenchmark, gSystem
20
21
# Create a new canvas, and customize it.
22
c1 =
TCanvas
(
"c1"
,
"Dynamic Filling Example"
, 200, 10, 700, 500)
23
c1.SetFillColor
(42)
24
frame =
c1.GetFrame
()
25
frame.SetFillColor
(21)
26
frame.SetBorderSize
(6)
27
frame.SetBorderMode
(-1)
28
29
# Create some histograms, a profile histogram and an ntuple
30
hpx =
TH1F
(
"hpx"
,
"This is the px distribution"
, 100, -4, 4)
31
hpxpy =
TH2F
(
"hpxpy"
,
"py vs px"
, 40, -4, 4, 40, -4, 4)
32
hprof =
TProfile
(
"hprof"
,
"Profile of pz versus px"
, 100, -4, 4, 0, 20)
33
ntuple =
TNtuple
(
"ntuple"
,
"Demo ntuple"
,
"px:py:pz:random:i"
)
34
35
# Set canvas/frame attributes.
36
hpx.SetFillColor
(48)
37
38
gBenchmark.Start
(
"hsimple"
)
39
40
# Fill histograms randomly.
41
for
i
in
range
(25000):
42
# Retrieve the generated values
43
px, py =
numpy.random.standard_normal
(size=2)
44
45
pz = px * px + py * py
46
random =
numpy.random.rand
()
47
48
# Fill histograms.
49
hpx.Fill
(px)
50
hpxpy.Fill
(px, py)
51
hprof.Fill
(px, pz)
52
ntuple.Fill
(px, py, pz, random, i)
53
54
# Update display every 1000 events.
55
if
i > 0
and
i % 1000 == 0:
56
hpx.Draw
()
57
58
c1.Modified
()
59
c1.Update
()
60
61
if
gSystem.ProcessEvents
():
# allow user interrupt
62
break
63
64
gBenchmark.Show
(
"hsimple"
)
65
66
# Create a new ROOT binary machine independent file.
67
# Note that this file may contain any kind of ROOT objects, histograms,
68
# pictures, graphics objects, detector geometries, tracks, events, etc.
69
with
TFile
(
"py-hsimple.root"
,
"RECREATE"
,
"Demo ROOT file with histograms"
)
as
hfile:
70
# Save all created objects in the file
71
hfile.WriteObject
(hpx)
72
hfile.WriteObject
(hpxpy)
73
hfile.WriteObject
(hprof)
74
# TNTuple is special because it is a TTree-derived class. To make sure all the
75
# dataset is properly written to disk, we connect it to the file and then
76
# we ask the ntuple to write all the information to the file itself.
77
ntuple.SetDirectory
(hfile)
78
ntuple.Write
()
79
80
c1.Modified
()
81
c1.Update
()
TRangeDynCast
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Definition
TCollection.h:358
ROOT::Detail::TRangeCast
Definition
TCollection.h:311
TCanvas
The Canvas class.
Definition
TCanvas.h:23
TFile
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition
TFile.h:131
TH1F
1-D histogram with a float per channel (see TH1 documentation)
Definition
TH1.h:879
TH2F
2-D histogram with a float per channel (see TH1 documentation)
Definition
TH2.h:307
TNtuple
A simple TTree restricted to a list of float variables only.
Definition
TNtuple.h:28
TProfile
Profile Histogram.
Definition
TProfile.h:32
tutorials
hsimple.py
ROOT master - Reference Guide Generated on Sun Sep 21 2025 15:08:07 (GVA Time) using Doxygen 1.10.0