← Back to branch summary

~madteam/mg5amcnlo/series2.0

« back to all changes in this revision

Viewing changes to madgraph/iolibs/template_files/matrix_madevent_v4.inc

  • Committer: olivier Mattelaer
  • Date: 2015-03-05 00:14:16 UTC
  • mfrom: (258.1.9 2.3)
  • mto: (258.8.1 2.3)
  • mto: This revision was merged to the branch mainline in revision 259.
  • Revision ID: olivier.mattelaer@uclouvain.be-20150305001416-y9mzeykfzwnl9t0j
partial merge

Show diffs side-by-side

added added

removed removed

Lines of Context:
madgraph/iolibs/template_files/matrix_madevent_v4.inc
10
10
C  
11
11
%(process_lines)s
12
12
C  
 
13
    use DiscreteSampler
13
14
    IMPLICIT NONE
14
15
C  
15
16
C CONSTANTS
30
31
C ARGUMENTS 
31
32
C  
32
33
    REAL*8 P(0:3,NEXTERNAL),ANS
 
34
c
 
35
c global (due to reading writting)
 
36
c 
 
37
    LOGICAL GOODHEL(NCOMB)
 
38
    INTEGER NTRY
 
39
    common/BLOCK_GOODHEL/NTRY,GOODHEL
33
40
C  
34
41
C LOCAL VARIABLES 
35
42
C  
36
 
    INTEGER NHEL(NEXTERNAL,NCOMB),NTRY
 
43
    INTEGER NHEL(NEXTERNAL,NCOMB)
37
44
    REAL*8 T,MATRIX%(proc_id)s
38
45
    REAL*8 R,SUMHEL,TS(NCOMB)
39
46
    INTEGER I,IDEN
40
47
    INTEGER IPROC,JC(NEXTERNAL),II
41
 
    LOGICAL GOODHEL(NCOMB)
42
48
    REAL*8 HWGT, XTOT, XTRY, XREJ, XR, YFRAC(0:NCOMB)
43
49
    INTEGER IDUM, NGOOD, IGOOD(NCOMB), JHEL, J, JJ
44
50
    REAL     XRAN1
59
65
    LOGICAL                    MULTI_CHANNEL
60
66
    COMMON/TO_MATRIX/ISUM_HEL, MULTI_CHANNEL
61
67
%(define_iconfigs_lines)s
62
 
    DATA NTRY,IDUM /0,-1/
 
68
    DATA IDUM /-1/
63
69
    DATA XTRY, XREJ, NGOOD /0,0,0/
64
70
    SAVE YFRAC, IGOOD, JHEL
65
 
    DATA GOODHEL/THEL*.FALSE./
66
71
%(helicity_lines)s
67
72
%(den_factor_line)s
68
73
C ----------
87
92
    DO I=1,NCOMB
88
93
       TS(I)=0d0
89
94
    ENDDO
90
 
    IF (ISUM_HEL .EQ. 0 .OR. NTRY .LE. MAXTRIES) THEN
 
95
 
 
96
!   If the helicity grid status is 0, this means that it is not yet initialized.
 
97
    IF (ISUM_HEL.EQ.0.or.(DS_get_dim_status('Helicity').eq.0)) THEN
91
98
        DO I=1,NCOMB
92
 
           IF (GOODHEL(I) .OR. NTRY .LE. MAXTRIES) THEN
 
99
           IF (GOODHEL(I) .OR. NTRY .LE. MAXTRIES.OR.(ISUM_HEL.NE.0)) THEN
93
100
               T=MATRIX%(proc_id)s(P ,NHEL(1,I),JC(1))            
94
101
             DO JJ=1,nincoming
95
102
               IF(POL(JJ).NE.1d0.AND.NHEL(JJ,I).EQ.INT(SIGN(1d0,POL(JJ)))) THEN
98
105
                 T=T*(2d0-ABS(POL(JJ)))
99
106
               ENDIF
100
107
             ENDDO
101
 
             ANS=ANS+T
 
108
                         IF (ISUM_HEL.NE.0) then
 
109
                           call DS_add_entry('Helicity',I,T)
 
110
                         endif
 
111
             ANS=ANS+DABS(T)
102
112
             TS(I)=T
103
113
           ENDIF
104
114
        ENDDO
105
 
        JHEL = 1
106
 
        IF(NTRY.LE.MAXTRIES)THEN
107
 
         DO I=1,NCOMB
108
 
            IF (.NOT.GOODHEL(I) .AND. (TS(I).GT.ANS*LIMHEL/NCOMB)) THEN
109
 
               GOODHEL(I)=.TRUE.
110
 
               NGOOD = NGOOD +1
111
 
               IGOOD(NGOOD) = I        
112
 
               print *,'Adding good helicity ',I,TS(I)/ANS
113
 
            ENDIF
114
 
         ENDDO
115
 
        ENDIF
116
 
        IF(NTRY.EQ.MAXTRIES)THEN
117
 
           ISUM_HEL=MIN(ISUM_HEL,NGOOD)
118
 
        ENDIF
 
115
        IF(NTRY.EQ.(MAXTRIES+1)) THEN
 
116
           call reset_cumulative_variable() ! avoid biais of the initialization
 
117
        ENDIF
 
118
        IF (ISUM_HEL.NE.0) then
 
119
!         We set HEL_PICKED to -1 here so that later on, the call to DS_add_point in dsample.f does not add anything to the grid since it was already done here.
 
120
                  HEL_PICKED = -1
 
121
!         For safety, hardset the helicity sampling jacobian to 0.0d0 to make sure it is not .
 
122
                  hel_jacobian   = 1.0d0
 
123
                  IF(DS_get_dim_status('Helicity').eq.1) then 
 
124
!           If we finished the initialization we can update the grid so as to start sampling over it.
 
125
!           However the grid will now be filled by dsample with different kind of weights (including pdf, flux, etc...) so by setting the grid_mode of the reference grid to 'initialization' we make sure it will be overwritten (as opposed to 'combined') by the running grid at the next update.
 
126
            CALL DS_UPDATE_GRID('Helicity')
 
127
                        CALL DS_SET_GRID_MODE('Helicity','init')
 
128
          endif
 
129
            ELSE
 
130
          JHEL = 1
 
131
          IF(NTRY.LE.MAXTRIES)THEN
 
132
           DO I=1,NCOMB
 
133
              IF (.NOT.GOODHEL(I) .AND. (TS(I).GT.ANS*LIMHEL/NCOMB)) THEN
 
134
                 GOODHEL(I)=.TRUE.
 
135
                 NGOOD = NGOOD +1
 
136
                 IGOOD(NGOOD) = I        
 
137
                 print *,'Adding good helicity ',I,TS(I)/ANS
 
138
              ENDIF
 
139
           ENDDO
 
140
          ENDIF
 
141
          IF(NTRY.EQ.MAXTRIES)THEN
 
142
             ISUM_HEL=MIN(ISUM_HEL,NGOOD)
 
143
          ENDIF
 
144
                endif
119
145
    ELSE              !RANDOM HELICITY
120
 
        DO J=1,ISUM_HEL
121
 
            JHEL=JHEL+1
122
 
            IF (JHEL .GT. NGOOD) JHEL=1
123
 
            HWGT = REAL(NGOOD)/REAL(ISUM_HEL)
124
 
            I = IGOOD(JHEL)
125
 
            T=MATRIX%(proc_id)s(P ,NHEL(1,I),JC(1))            
 
146
       
 
147
C           The helicity configuration was chosen already by genps and put in a common block defined in genps.inc.
 
148
            I = HEL_PICKED
 
149
                        
 
150
                        T=MATRIX%(proc_id)s(P ,NHEL(1,I),JC(1))            
126
151
            DO JJ=1,nincoming
127
152
              IF(POL(JJ).NE.1d0.AND.NHEL(JJ,I).EQ.INT(SIGN(1d0,POL(JJ)))) THEN
128
153
                T=T*ABS(POL(JJ))
130
155
                T=T*(2d0-ABS(POL(JJ)))
131
156
              ENDIF
132
157
            ENDDO
133
 
            ANS=ANS+T*HWGT
134
 
            TS(I)=T*HWGT
135
 
        ENDDO
136
 
        IF (ISUM_HEL .EQ. 1) THEN
137
 
            WRITE(HEL_BUFF,'(20i5)')(NHEL(II,I),II=1,NEXTERNAL)
138
 
        ENDIF
 
158
c           Always one helicity at a time
 
159
            ANS = T
 
160
c           Include the Jacobian from helicity sampling
 
161
            ANS = ANS * hel_jacobian
 
162
            WRITE(HEL_BUFF,'(20i5)')(NHEL(II,I),II=1,NEXTERNAL)  
139
163
    ENDIF
140
 
    IF (ISUM_HEL .NE. 1) THEN
 
164
    IF (ISUM_HEL .NE. 1.or.(HEL_PICKED.eq.-1)) THEN
141
165
    R=XRAN1(IDUM)*ANS
142
166
    SUMHEL=0d0
143
167
    DO I=1,NCOMB
144
168
       SUMHEL=SUMHEL+TS(I)
145
169
       IF(R.LT.SUMHEL)THEN
146
170
          WRITE(HEL_BUFF,'(20i5)')(NHEL(II,I),II=1,NEXTERNAL)
 
171
          ANS=DSIGN(ANS,TS(I))            
147
172
          GOTO 10
148
173
       ENDIF
149
174
    ENDDO