← Back to branch summary

~maddevelopers/mg5amcnlo/3.0.1

« back to all changes in this revision

Viewing changes to Template/NLO/FixedOrderAnalysis/analysis_root_template.f

  • Committer: Marco Zaro
  • Date: 2014-01-27 16:54:10 UTC
  • mfrom: (78.124.55 MG5_aMC_2.1)
  • Revision ID: marco.zaro@gmail.com-20140127165410-5lma8c2hzbzm426j
merged with lp:~maddevelopers/madgraph5/MG5_aMC_2.1 r 267

Show diffs side-by-side

added added

removed removed

Lines of Context:
Template/NLO/FixedOrderAnalysis/analysis_root_template.f
 
1
c
 
2
c This file contains the default histograms for fixed order runs: it
 
3
c only plots the total rate as an example. It can be used as a template
 
4
c to make distributions for other observables.
 
5
c
 
6
c It generates histograms in root format, by using fortran front-end
 
7
c interfaces originally written for MC@NLO by Wouter Verkerke
 
8
c (see the rbook* files)
 
9
c
 
10
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
11
      subroutine analysis_begin(nwgt,weights_info)
 
12
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
13
c This subroutine is called once at the start of each run. Here the
 
14
c histograms should be declared. 
 
15
c
 
16
c Declare the histograms using 'rbook'.
 
17
c     o) The first argument is an integer that labels the histogram. In
 
18
c     the analysis_end and analysis_fill subroutines this label is used
 
19
c     to keep track of the histogram. The label should be a number
 
20
c     between 1 and NPLOTS/4=1250 (can be increased in dbook.inc).
 
21
c     o) The second argument is a string that will apear above the
 
22
c     histogram.
 
23
c     o) The third, forth and fifth arguments are the bin size, the
 
24
c     lower edge of the first bin and the upper edge of the last
 
25
c     bin. 
 
26
c     o) When including scale and/or PDF uncertainties, declare a
 
27
c     histogram for each weight, and compute the uncertainties from the
 
28
c     final set of histograms
 
29
c
 
30
      implicit none
 
31
c When including scale and/or PDF uncertainties the total number of
 
32
c weights considered is nwgt
 
33
      integer nwgt
 
34
c In the weights_info, there is an text string that explains what each
 
35
c weight will mean. The size of this array of strings is equal to nwgt.
 
36
      character*(*) weights_info(*)
 
37
c Local variables
 
38
      integer kk,l,nwgt_analysis
 
39
      common/c_analysis/nwgt_analysis
 
40
c Initialize the histogramming package (rbook): 
 
41
      call open_root_file()
 
42
c Fill the c_analysis common block with the number of weights that will
 
43
c be computed
 
44
      nwgt_analysis=nwgt
 
45
c
 
46
c loop over all the weights that are computed (depends on run_card
 
47
c parameters do_rwgt_scale and do_rwgt_pdf):
 
48
      do kk=1,nwgt_analysis
 
49
c make sure that there is a separate histogram initialized for each
 
50
c weight
 
51
         l=(kk-1)*2
 
52
c declare (i.e. book) the histograms
 
53
         call rbook(l+1,'total rate      '//weights_info(kk),
 
54
     &        1.0d0,0.5d0,5.5d0)
 
55
         call rbook(l+2,'total rate Born '//weights_info(kk),
 
56
     &        1.0d0,0.5d0,5.5d0)
 
57
      enddo
 
58
      return
 
59
      end
 
60
 
 
61
 
 
62
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
63
      subroutine analysis_end(xnorm)
 
64
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
65
c This subroutine is called once at the end of the run. Here the
 
66
c histograms are written to disk. Note that this is done for each
 
67
c integration channel separately. There is an external script that will
 
68
c read the root files in each of the integration channels and
 
69
c combines them by summing all the bins in a final single root
 
70
c file to be put in the Events/run_XX directory.
 
71
c      
 
72
c
 
73
c The user is NOT SUPPOSED TO MODIFY THIS ROUTINE, except for making
 
74
c sure that the calls to "ropera" below are relevant to all
 
75
c histograms previously booked
 
76
c
 
77
      implicit none
 
78
      double precision xnorm
 
79
c Local variables
 
80
      integer kk,l,nwgt_analysis
 
81
      common/c_analysis/nwgt_analysis
 
82
c Do not touch the following lines. These lines make sure that the
 
83
c histograms will have the correct overall normalisation: cross section
 
84
c (in pb) per bin.
 
85
      do kk=1,nwgt_analysis
 
86
         l=(kk-1)*2
 
87
         call ropera(l+1,'+',l+1,l+1,xnorm,0.d0)
 
88
         call ropera(l+2,'+',l+2,l+2,xnorm,0.d0)
 
89
      enddo
 
90
      call close_root_file
 
91
      return
 
92
      end
 
93
 
 
94
 
 
95
 
 
96
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
97
      subroutine analysis_fill(p,istatus,ipdg,wgts,ibody)
 
98
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 
99
c This subroutine is called for each n-body and (n+1)-body configuration
 
100
c that passes the generation cuts. Here the histrograms are filled.
 
101
      implicit none
 
102
c This includes the 'nexternal' parameter that labels the number of
 
103
c particles in the (n+1)-body process
 
104
      include 'nexternal.inc'
 
105
c This is an array which is '-1' for initial state and '1' for final
 
106
c state particles
 
107
      integer istatus(nexternal)
 
108
c This is an array with (simplified) PDG codes for the particles. Note
 
109
c that channels that are combined (i.e. they have the same matrix
 
110
c elements) are given only 1 set of PDG codes. This means, e.g., that
 
111
c when using a 5-flavour scheme calculation (massless b quark), no
 
112
c b-tagging can be applied.
 
113
      integer iPDG(nexternal)
 
114
c The array of the momenta and masses of the initial and final state
 
115
c particles in the lab frame. The format is "E, px, py, pz, mass", while
 
116
c the second dimension loops over the particles in the process. Note
 
117
c that these are the (n+1)-body particles; for the n-body there is one
 
118
c momenta equal to all zero's (this is not necessarily the last particle
 
119
c in the list). If one uses IR-safe obserables only, there should be no
 
120
c difficulty in using this.
 
121
      double precision p(0:4,nexternal)
 
122
c The weight of the current phase-space point is wgts(1). If scale
 
123
c and/or PDF uncertainties are included through reweighting, the rest of
 
124
c the array contains the list of weights in the same order as described
 
125
c by the weigths_info strings in analysis_begin
 
126
      double precision wgts(*)
 
127
c The ibody variable is:
 
128
c     ibody=1 : (n+1)-body contribution
 
129
c     ibody=2 : n-body contribution (excluding the Born)
 
130
c     ibody=3 : Born contribution
 
131
c The histograms need to be filled for all these contribution to get the
 
132
c physics NLO results. (Note that the adaptive phase-space integration
 
133
c is optimized using the sum of the contributions, therefore plotting
 
134
c them separately might lead to larger than expected statistical
 
135
c fluctuations).
 
136
      integer ibody
 
137
c local variables
 
138
      double precision wgt,var
 
139
      integer kk,l,nwgt_analysis
 
140
      common/c_analysis/nwgt_analysis
 
141
c
 
142
c Fill the histograms here using a call to the rfill() subroutine. The
 
143
c first argument is the histogram label, the second is the numerical
 
144
c value of the variable to plot for the current phase-space point and
 
145
c the final argument is the weight of the current phase-space point.
 
146
      var=1d0
 
147
c loop over all the weights that are computed (depends on run_card
 
148
c parameters do_rwgt_scale and do_rwgt_pdf):
 
149
      do kk=1,nwgt_analysis
 
150
         wgt=wgts(kk)
 
151
         l=(kk-1)*2
 
152
c always fill the total rate
 
153
         call rfill(l+1,var,wgt)
 
154
c only fill the total rate for the Born when ibody=3
 
155
         if (ibody.eq.3) call rfill(l+2,var,wgt)
 
156
      enddo
 
157
      return
 
158
      end