ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
zdemo.py
Go to the documentation of this file.
1 # This macro is an example of graphs in log scales with annotations.
2 #
3 # The begin_html <a href="gif/zdemo.gif" >presented results</a> end_html
4 # are predictions of invariant cross-section of Direct Photons produced
5 # at RHIC energies, based on the universality of scaling function H(z).
6 #
7 # Authors: Michael Tokarev and Elena Potrebenikova (JINR Dubna)
8 #
9 # These Figures were published in JINR preprint E2-98-64, Dubna,
10 # 1998 and submitted to CPC.
11 #
12 # Note that the way greek symbols, super/subscripts are obtained
13 # illustrate the current limitations of Root in this area.
14 #
15 
16 import ROOT
17 from array import array
18 from math import *
19 
20 NMAX = 20
21 Z = array( 'f', [0.]*NMAX )
22 HZ = array( 'f', [0.]*NMAX )
23 PT = array( 'f', [0.]*NMAX )
24 INVSIG = array( 'f', [0.]*NMAX )
25 
26 NLOOP = 0
27 saves = {}
28 
29 #_______________________________________________________________________________
30 def zdemo():
31  global NLOOP
32  global Z, HZ, PT, INVSIG
33  global saves
34 
35  # Create a new canvas.
36  c1 = ROOT.TCanvas( 'zdemo', 'Monte Carlo Study of Z scaling', 10, 40, 800, 600 )
37  c1.Range( 0, 0, 25, 18 )
38  c1.SetFillColor( 40 )
39  saves[ 'c1' ] = c1 # prevent deteletion at end of zdemo
40 
41  pl = ROOT.TPaveLabel( 1, 16.3, 24, 17.5,
42  'Z-scaling of Direct Photon Productions in pp Collisions at RHIC Energies', 'br' )
43  pl.SetFillColor(18)
44  pl.SetTextFont(32)
45  pl.SetTextColor(49)
46  pl.Draw()
47  saves[ 'pl' ] = pl
48 
49  t = ROOT.TLatex()
50  t.SetTextFont(32)
51  t.SetTextColor(1)
52  t.SetTextSize(0.03)
53  t.SetTextAlign(12)
54  t.DrawLatex( 3.1, 15.5, 'M.Tokarev, E.Potrebenikova ')
55  t.DrawLatex( 14., 15.5, 'JINR preprint E2-98-64, Dubna, 1998 ')
56  saves[ 't' ] = t
57 
58  pad1 = ROOT.TPad( 'pad1', 'This is pad1', 0.02, 0.02, 0.48, 0.83, 33 )
59  pad2 = ROOT.TPad( 'pad2', 'This is pad2', 0.52, 0.02, 0.98, 0.83, 33 )
60 
61  pad1.Draw()
62  pad2.Draw()
63 
64  saves[ 'pad1' ] = pad1; saves[ 'pad2' ] = pad2
65 
66 #
67 # Cross-section of direct photon production in pp collisions at 500 GeV vs Pt
68 #
69  energ = 63
70  dens = 1.766
71  tgrad = 90.
72  ptmin = 4.
73  ptmax = 24.
74  delp = 2.
75  hz_calc( energ, dens, tgrad, ptmin, ptmax, delp )
76  pad1.cd()
77  pad1.Range( -0.255174, -19.25, 2.29657, -6.75 )
78  pad1.SetLogx()
79  pad1.SetLogy()
80 
81  # create a 2-d histogram to define the range
82  pad1.DrawFrame( 1, 1e-18, 110, 1e-8 )
83  pad1.GetFrame().SetFillColor( 19 )
84  t = ROOT.TLatex()
85  t.SetNDC()
86  t.SetTextFont( 62 )
87  t.SetTextColor( 36 )
88  t.SetTextSize( 0.08 )
89  t.SetTextAlign( 12 )
90  t.DrawLatex( 0.6, 0.85, 'p - p' )
91 
92  t.SetTextSize( 0.05 )
93  t.DrawLatex( 0.6, 0.79, 'Direct #gamma' )
94  t.DrawLatex( 0.6, 0.75, '#theta = 90^{o}' )
95 
96  t.DrawLatex( 0.20, 0.45, 'Ed^{3}#sigma/dq^{3}' )
97  t.DrawLatex( 0.18, 0.40, '(barn/Gev^{2})' )
98 
99  t.SetTextSize( 0.045 )
100  t.SetTextColor( ROOT.kBlue )
101  t.DrawLatex( 0.22, 0.260, '#sqrt{s} = 63(GeV)' )
102  t.SetTextColor( ROOT.kRed )
103  t.DrawLatex( 0.22, 0.205,'#sqrt{s} = 200(GeV)' )
104  t.SetTextColor( 6 )
105  t.DrawLatex( 0.22, 0.15, '#sqrt{s} = 500(GeV)' )
106 
107  t.SetTextSize( 0.05 )
108  t.SetTextColor( 1 )
109  t.DrawLatex( 0.6, 0.06, 'q_{T} (Gev/c)' )
110  saves[ 't2' ] = t # note the label that is used!
111 
112  gr1 = ROOT.TGraph( NLOOP, PT, INVSIG )
113 
114  gr1.SetLineColor( 38 )
115  gr1.SetMarkerColor( ROOT.kBlue )
116  gr1.SetMarkerStyle( 21 )
117  gr1.SetMarkerSize( 1.1 )
118  gr1.Draw( 'LP' )
119  saves[ 'gr1' ] = gr1
120 
121 #
122 # Cross-section of direct photon production in pp collisions at 200 GeV vs Pt
123 #
124 
125  energ = 200
126  dens = 2.25
127  tgrad = 90.
128  ptmin = 4.
129  ptmax = 64.
130  delp = 6.
131  hz_calc( energ, dens, tgrad, ptmin, ptmax, delp )
132 
133  gr2 = ROOT.TGraph( NLOOP, PT, INVSIG )
134  gr2.SetLineColor( 38 )
135  gr2.SetMarkerColor( ROOT.kRed )
136  gr2.SetMarkerStyle( 29 )
137  gr2.SetMarkerSize( 1.5 )
138  gr2.Draw( 'LP' )
139  saves[ 'gr2' ] = gr2
140 
141 #
142 # Cross-section of direct photon production in pp collisions at 500 GeV vs Pt
143 #
144  energ = 500
145  dens = 2.73
146  tgrad = 90.
147  ptmin = 4.
148  ptmax = 104.
149  delp = 10.
150  hz_calc( energ, dens, tgrad, ptmin, ptmax, delp )
151 
152  gr3 = ROOT.TGraph( NLOOP, PT, INVSIG )
153 
154  gr3.SetLineColor( 38 )
155  gr3.SetMarkerColor( 6 )
156  gr3.SetMarkerStyle( 8 )
157  gr3.SetMarkerSize( 1.1 )
158  gr3.Draw( 'LP' )
159  saves[ 'gr3' ] = gr3
160 
161  dum = array( 'f', [0.] )
162  graph = ROOT.TGraph( 1, dum, dum )
163  graph.SetMarkerColor( ROOT.kBlue )
164  graph.SetMarkerStyle( 21 )
165  graph.SetMarkerSize( 1.1 )
166  graph.SetPoint( 0, 1.7, 1.e-16 )
167  graph.Draw( 'LP' )
168  saves[ 'graph' ] = graph
169 
170  graph = ROOT.TGraph( 1, dum, dum )
171  graph.SetMarkerColor( ROOT.kRed )
172  graph.SetMarkerStyle( 29 )
173  graph.SetMarkerSize( 1.5 )
174  graph.SetPoint( 0, 1.7, 2.e-17 )
175  graph.Draw( 'LP' )
176  saves[ 'graph2' ] = graph # note the label that is used!
177 
178  graph = ROOT.TGraph( 1, dum, dum )
179  graph.SetMarkerColor( 6 )
180  graph.SetMarkerStyle( 8 )
181  graph.SetMarkerSize( 1.1 )
182  graph.SetPoint( 0, 1.7, 4.e-18)
183  graph.Draw( 'LP' )
184  saves[ 'graph3' ] = graph # note the label that is used!
185 
186  pad2.cd()
187  pad2.Range( -0.43642, -23.75, 3.92778, -6.25 )
188  pad2.SetLogx()
189  pad2.SetLogy()
190 
191  pad2.DrawFrame( 1, 1e-22, 3100, 1e-8 )
192  pad2.GetFrame().SetFillColor( 19 )
193 
194  gr = ROOT.TGraph( NLOOP, Z, HZ )
195  gr.SetTitle( 'HZ vs Z' )
196  gr.SetFillColor( 19 )
197  gr.SetLineColor( 9 )
198  gr.SetMarkerColor( 50 )
199  gr.SetMarkerStyle( 29 )
200  gr.SetMarkerSize( 1.5 )
201  gr.Draw( 'LP' )
202  saves[ 'gr' ] = gr
203 
204  t = ROOT.TLatex()
205  t.SetNDC()
206  t.SetTextFont( 62 )
207  t.SetTextColor( 36 )
208  t.SetTextSize( 0.08 )
209  t.SetTextAlign( 12 )
210  t.DrawLatex( 0.6, 0.85, 'p - p' )
211 
212  t.SetTextSize( 0.05 )
213  t.DrawLatex( 0.6, 0.79, 'Direct #gamma' )
214  t.DrawLatex( 0.6, 0.75, '#theta = 90^{o}' )
215 
216  t.DrawLatex( 0.70, 0.55, 'H(z)' )
217  t.DrawLatex( 0.68, 0.50, '(barn)' )
218 
219  t.SetTextSize( 0.045 )
220  t.SetTextColor( 46 )
221  t.DrawLatex( 0.20, 0.30, '#sqrt{s}, GeV' )
222  t.DrawLatex( 0.22, 0.26, '63' )
223  t.DrawLatex( 0.22, 0.22, '200' )
224  t.DrawLatex( 0.22, 0.18, '500' )
225 
226  t.SetTextSize( 0.05 )
227  t.SetTextColor( 1 )
228  t.DrawLatex( 0.88, 0.06, 'z' )
229  saves[ 't3' ] = t # note the label that is used!
230 
231  c1.Modified()
232  c1.Update()
233 
234 #_______________________________________________________________________________
235 def hz_calc( ENERG, DENS, TGRAD, PTMIN, PTMAX, DELP ):
236  global NLOOP
237  global Z, HZ, PT, INVSIG
238 
239  CSEFT= 1.
240  GM1 = 0.00001
241  GM2 = 0.00001
242  A1 = 1.
243  A2 = 1.
244  ALX = 2.
245  BETA = 1.
246  KF1 = 8.E-7
247  KF2 = 5.215
248 
249  MN = 0.9383
250  DEGRAD=0.01745329
251 
252  # print 'ENR= %f DENS= %f PTMIN= %f PTMAX= %f DELP= %f ' % (ENERG,DENS,PTMIN,PTMAX,DELP)
253 
254  DNDETA= DENS
255  MB1 = MN*A1
256  MB2 = MN*A2
257  EB1 = ENERG/2.*A1
258  EB2 = ENERG/2.*A2
259  M1 = GM1
260  M2 = GM2
261  THET = TGRAD*DEGRAD
262  NLOOP = int((PTMAX-PTMIN)/DELP)
263 
264  for I in range(NLOOP):
265  PT[I]=PTMIN+I*DELP
266  PTOT = PT[I]/sin(THET)
267 
268  ETOT = sqrt(M1*M1 + PTOT*PTOT)
269  PB1 = sqrt(EB1*EB1 - MB1*MB1)
270  PB2 = sqrt(EB2*EB2 - MB2*MB2)
271  P2P3 = EB2*ETOT+PB2*PTOT*cos(THET)
272  P1P2 = EB2*EB1+PB2*PB1
273  P1P3 = EB1*ETOT-PB1*PTOT*cos(THET)
274 
275  X1 = P2P3/P1P2
276  X2 = P1P3/P1P2
277  Y1 = X1+sqrt(X1*X2*(1.-X1)/(1.-X2))
278  Y2 = X2+sqrt(X1*X2*(1.-X2)/(1.-X1))
279 
280  S = (MB1*MB1)+2.*P1P2+(MB2*MB2)
281  SMIN = 4.*((MB1*MB1)*(X1*X1) +2.*X1*X2*P1P2+(MB2*MB2)*(X2*X2))
282  SX1 = 4.*( 2*(MB1*MB1)*X1+2*X2*P1P2)
283  SX2 = 4.*( 2*(MB2*MB2)*X2+2*X1*P1P2)
284  SX1X2= 4.*(2*P1P2)
285  DELM = pow((1.-Y1)*(1.-Y2),ALX)
286 
287  Z[I] = sqrt(SMIN)/DELM/pow(DNDETA,BETA)
288 
289  Y1X1 = 1. +X2*(1-2.*X1)/(2.*(Y1-X1)*(1.-X2))
290  Y1X2 = X1*(1-X1)/(2.*(Y1-X1)*(1.-X2)*(1.-X2))
291  Y2X1 = X2*(1-X2)/(2.*(Y2-X2)*(1.-X1)*(1.-X1))
292  Y2X2 = 1. +X1*(1-2.*X2)/(2.*(Y2-X2)*(1.-X1))
293  Y2X1X2= Y2X1*( (1.-2.*X2)/(X2*(1-X2)) -( Y2X2-1.)/(Y2-X2))
294  Y1X1X2= Y1X2*( (1.-2.*X1)/(X1*(1-X1)) -( Y1X1-1.)/(Y1-X1))
295 
296  KX1=-DELM*(Y1X1*ALX/(1.-Y1) + Y2X1*ALX/(1.-Y2))
297  KX2=-DELM*(Y2X2*ALX/(1.-Y2) + Y1X2*ALX/(1.-Y1))
298  ZX1=Z[I]*(SX1/(2.*SMIN)-KX1/DELM)
299  ZX2=Z[I]*(SX2/(2.*SMIN)-KX2/DELM)
300 
301  H1=ZX1*ZX2
302 
303  HZ[I]=KF1/pow(Z[I],KF2)
304  INVSIG[I]=(HZ[I]*H1*16.)/S
305 
306 
307 # run if loaded as script
308 if __name__ == '__main__':
309  zdemo()
zdemo.zdemo
def zdemo
Definition: zdemo.py:30
cos
double cos(double)
zdemo
Definition: zdemo.py:1
sqrt
double sqrt(double)
zdemo.hz_calc
def hz_calc
Definition: zdemo.py:235
pow
double pow(double, double)
sin
double sin(double)
SetFillColor
h1 SetFillColor(kGreen)