EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
poldsigma.F
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file poldsigma.F
1 *72*********************************************************************
2 
3  FUNCTION poldsigma(XP)
4 
5  IMPLICIT NONE
6 C...Differential cross section for first order QCD processes.
7 
8 
9 *
10 * to avoid variable conflictions, a second keep element is necessary
11 * with the same common block name (see LPTOU2)
12 *
13 *
14 * to avoid variable conflictions, a second keep element is necessary
15 * with the same common block name (see LPTOU2)
16 *
17  COMMON /leptou/ cut(14),lst(40),parl(30),
18  & x,y,w2,q2,u
19  REAL cut,parl,x,y,w2,q2,u
20  INTEGER lst
21  SAVE /leptou/
22 
23  COMMON /linter/ pari(50),ewqc(2,2,8),qc(8),zl(2,4),zq(2,8),pq(17)
24  REAL pari,ewqc,qc,zl,zq,pq
25  SAVE /linter/
26 
27  COMMON /lintrl/ psave(3,4,5),ksave(4),xmin,xmax,ymin,ymax,
28  &q2min,q2max,w2min,w2max,ilep,inu,ig,iz
29  REAL psave,xmin,xmax,ymin,ymax,q2min,q2max,w2min,w2max
30  INTEGER ksave,ilep,inu,ig,iz
31  SAVE /lintrl/
32 
33  INTEGER nlupdm,nplbuf
34  parameter(nlupdm=4000,nplbuf=5)
35  common/lujets/n,k(nlupdm,5),p(nlupdm,nplbuf),v(nlupdm,5)
36  INTEGER n,k
37  REAL p,v
38  SAVE /lujets/
39 
40  common/ludat1/mstu(200),paru(200),mstj(200),parj(200)
41  INTEGER mstu,mstj
42  REAL paru,parj
43  SAVE /ludat1/
44 
45 
46  REAL poldsigma,xp,dqcdi,poldqcdi,ulmass
47  INTEGER i,ih,il,iu,ip
48  REAL pqh,amu,xi,zpmin,zpmax,xpq,wq,wqb,tq,tqb,t1,s13,sgn,xdpq,sig
49 *...Added array XDPQ to store delta parton distributions.
50  dimension xpq(-6:6),xdpq(-6:6),pqh(17,2)
51 *___
52  poldsigma=0.
53  DO 10 i=1,17
54  pqh(i,1)=0.
55  pqh(i,2)=0.
56  10 pq(i)=0.
57 
58  mstj(93)=1
59  amu=ulmass(2)
60 *PEPSI>>
61  il=1
62  iu=3
63 *PEPSI<<
64  xi=x/xp
65 C...Scheme for ME cutoff: W2, Q2, mixed
66  IF(lst(20).LE.1) THEN
67  zpmin=(1.-x)*xp/(xp-x)*parl(27)
68  ELSEIF(lst(20).EQ.2) THEN
69  zpmin=x*xp/(xp-x)*parl(27)
70  ELSEIF(lst(20).GE.3.AND.lst(20).LE.5) THEN
71  zpmin=parl(27)
72  ELSEIF(lst(20).GE.6) THEN
73  zpmin=parl(8)
74  ENDIF
75  IF(zpmin.LE.0..OR.zpmin.GE.0.5) RETURN
76  zpmax=1.d0-dble(zpmin)
77  CALL lnstrf(xi,q2,xpq)
78  CALL dnstrf(xi,q2,xdpq)
79  IF(lst(24).EQ.3) goto 3000
80 
81 C...Gluon bremsstrahlung process, i.e. qg-event.
82  2000 DO 2500 ip=il,iu
83  sig=dqcdi(1,ip,xp,zpmin,zpmax)
84 *PEPSI>>
85  IF (ip.EQ.3) THEN
86  sig=poldqcdi(1,ip,xp,zpmin,zpmax)
87  ENDIF
88 *PEPSI<<
89  sgn=sign(1.,5.-2.*ip)
90  DO 2300 ih=1,2
91  IF(ih.EQ.1) THEN
92  IF(parl(6).GT.0.99) goto 2300
93  IF(lst(32).EQ.0.AND.lst(30).NE.-1) goto 2300
94  ELSEIF(ih.EQ.2) THEN
95  IF(parl(6).LT.-0.99) goto 2300
96  IF(lst(32).EQ.0.AND.lst(30).NE.1) goto 2300
97  ENDIF
98  IF(lst(32).NE.0) lst(30)=sign(1.,ih-1.5)
99  IF(lst(23).NE.2) THEN
100  DO 2100 i=1,lst(12)
101  wq=xpq(i)*sig*(ewqc(1,ih,i)+sgn*ewqc(2,ih,i))
102  wqb=xpq(-i)*sig*sgn*(ewqc(1,ih,i)+sgn*ewqc(2,ih,i))
103 *...EM polarized case: if IP=3 use proper formula for
104 *...the interference part of cross section.
105 *PEPSI>>
106  IF(ip.EQ.3)THEN
107  wq=lst(40)*xdpq(i)*sig*(ewqc(1,ih,i)+ewqc(2,ih,i))
108  wqb=lst(40)*xdpq(-i)*sig*(ewqc(1,ih,i)+ewqc(2,ih,i))
109  ENDIF
110 *PEPSI<<
111 C...Include y-dependence.
112  wq=wq*pari(23+ip)
113  wqb=wqb*pari(23+ip)
114  pqh(i,ih)=pqh(i,ih)+wq
115  pqh(i+lst(12),ih)=pqh(i+lst(12),ih)+wqb
116  pqh(17,ih)=pqh(17,ih)+wq+wqb
117  2100 CONTINUE
118  ELSEIF(lst(23).EQ.2) THEN
119 C...Zero CC cross-section for one helicity state.
120  IF(ksave(1).LT.0.AND.ih.EQ.1
121  & .OR.ksave(1).GT.0.AND.ih.EQ.2) goto 2300
122  IF(ip.EQ.3) THEN
123  tq=-lst(30)
124  tqb=-tq
125  ELSE
126  tq=1.
127  tqb=1.
128  ENDIF
129  DO 2200 i=1,lst(12)
130  t1=-k(3,2)*qc(i)
131  IF(t1.GT.0) THEN
132  wq=xpq(i)*sig*tq
133  wqb=0.
134  ELSE
135  wqb=xpq(-i)*sig*tqb
136  wq=0.
137  ENDIF
138 C...Include y-dependence.
139  wq=wq*pari(23+ip)
140  wqb=wqb*pari(23+ip)
141  pqh(i,ih)=pqh(i,ih)+wq
142  pqh(i+lst(12),ih)=pqh(i+lst(12),ih)+wqb
143  pqh(17,ih)=pqh(17,ih)+wq+wqb
144  2200 CONTINUE
145  ENDIF
146  2300 CONTINUE
147  2500 CONTINUE
148  DO 2600 i=1,17
149  2600 pq(i)=(1.-parl(6))/2.*pqh(i,1)+(1.+parl(6))/2.*pqh(i,2)
150  ih=1
151  IF(lst(30).EQ.1) ih=2
152  IF(lst(32).EQ.0) THEN
153 C...Simulation: cross section for chosen helicity state only.
154  poldsigma=pqh(17,ih)
155  ELSEIF(lst(19).EQ.-1) THEN
156 C...Integration event-by-event: normalize and include alpha_s*1/(1-xp)
157  poldsigma=pqh(17,ih)/pari(20)*parl(25)/(1.-xp)
158 C...Max of dsigma/dxp for L- and R-handed lepton.
159  IF(pqh(17,1).GT.pari(15)) pari(15)=pqh(17,1)
160  IF(pqh(17,2).GT.pari(16)) pari(16)=pqh(17,2)
161  ELSE
162 C...Integration: normalize and include alpha_s*1/(1-xp)
163  poldsigma=pqh(17,ih)/pari(20)*parl(25)/(1.-xp)
164  IF(lst(17).EQ.0.AND.lst(40).EQ.0) THEN
165 C...Fixed beam energy, max of dsigma/dxp for L- and R-handed lepton.
166  IF(pqh(17,1).GT.pari(15)) pari(15)=pqh(17,1)
167  IF(pqh(17,2).GT.pari(16)) pari(16)=pqh(17,2)
168  ELSE
169 C...Variable beam energy or polarizzed case, max of dsigma/dxp
170 *...for S, T, I contributions.
171  IF(pq(17)/pari(23+lst(32)).GT.pari(14+lst(32)))
172  & pari(14+lst(32))=pq(17)/pari(23+lst(32))
173  ENDIF
174  ENDIF
175  RETURN
176 C...Boson-gluon fusion, i.e. qq-event.
177  3000 s13=q2*(1.-xp)/xp
178  IF(s13.LT.4.*amu**2) RETURN
179  DO 3500 ip=il,iu
180 *PEPSI>>
181 *HI SIG=XPQ(0)*DQCDI(2,IP,XP,ZPMIN,ZPMAX)
182  IF(ip.EQ.3.) THEN
183  sig=float(lst(40))*xdpq(0)*poldqcdi(2,ip,xp,zpmin,zpmax)
184  ELSE
185  sig=xpq(0)*dqcdi(2,ip,xp,zpmin,zpmax)
186  ENDIF
187 *PEPSI<<
188  DO 3300 ih=1,2
189  IF(ih.EQ.1) THEN
190  IF(parl(6).GT.0.99) goto 3300
191  IF(lst(32).EQ.0.AND.lst(30).NE.-1) goto 3300
192  ELSEIF(ih.EQ.2) THEN
193  IF(parl(6).LT.-0.99) goto 3300
194  IF(lst(32).EQ.0.AND.lst(30).NE.1) goto 3300
195  ENDIF
196  IF(lst(32).NE.0) lst(30)=sign(1.,ih-1.5)
197  IF(lst(23).NE.2) THEN
198  DO 3100 i=1,lst(13)
199  mstj(93)=1
200  IF(s13.LT.4.*ulmass(i)**2) goto 3100
201  wq=sig/2.*(ewqc(1,ih,i)+ewqc(2,ih,i))
202  wqb=wq
203 C...Include y-dependence.
204  wq=wq*pari(23+ip)
205  wqb=wqb*pari(23+ip)
206  pqh(i,ih)=pqh(i,ih)+wq
207  pqh(i+lst(13),ih)=pqh(i+lst(13),ih)+wqb
208  pqh(17,ih)=pqh(17,ih)+wq+wqb
209  3100 CONTINUE
210  ELSEIF(lst(23).EQ.2) THEN
211 C...Zero CC cross-section for one helicity state.
212  IF(ksave(1).LT.0.AND.ih.EQ.1
213  & .OR.ksave(1).GT.0.AND.ih.EQ.2) goto 3300
214  DO 3200 i=1,lst(13)
215  mstj(93)=1
216  IF(s13.LT.(amu+ulmass(i))**2) goto 3200
217  IF(k(3,2)*qc(i).LT.0) THEN
218  wq=sig
219  wqb=0.
220  ELSE
221  wqb=sig
222  wq=0.
223  ENDIF
224 C...Include y-dependence.
225  wq=wq*pari(23+ip)
226  wqb=wqb*pari(23+ip)
227  pqh(i,ih)=pqh(i,ih)+wq
228  pqh(i+lst(13),ih)=pqh(i+lst(13),ih)+wqb
229  pqh(17,ih)=pqh(17,ih)+wq+wqb
230  3200 CONTINUE
231  ENDIF
232  3300 CONTINUE
233  3500 CONTINUE
234  DO 3600 i=1,17
235  3600 pq(i)=(1.-parl(6))/2.*pqh(i,1)+(1.+parl(6))/2.*pqh(i,2)
236  ih=1
237  IF(lst(30).EQ.1) ih=2
238  IF(lst(32).EQ.0) THEN
239 C...Simulation: cross section for chosen helicity state only.
240  poldsigma=pqh(17,ih)
241  ELSEIF(lst(19).EQ.-1) THEN
242 C...Integration event-by-event: normalize and include alpha_s*1/(1-xp)
243  poldsigma=pqh(17,ih)/pari(20)*parl(25)/(1.-xp)
244 C...Max of dsigma/dxp for L- and R-handed lepton.
245  IF(pqh(17,1).GT.pari(18)) pari(18)=pqh(17,1)
246  IF(pqh(17,2).GT.pari(19)) pari(19)=pqh(17,2)
247  ELSE
248 C...Integration for grid: normalize and include alpha_s*1/(1-xp)
249  poldsigma=pq(17)/pari(20)*parl(25)/(1.-xp)
250  IF(lst(17).EQ.0.AND.lst(40).EQ.0) THEN
251 *---
252 C...Fixed beam energy, max of dsigma/dxp for L- and R-handed lepton.
253  IF(pqh(17,1).GT.pari(18)) pari(18)=pqh(17,1)
254  IF(pqh(17,2).GT.pari(19)) pari(19)=pqh(17,2)
255  ELSE
256 C...Variable beam energy, max of dsigma/dxp for S, T, I contributions.
257  IF(pq(17)/pari(23+lst(32)).GT.pari(17+lst(32)))
258  & pari(17+lst(32))=pq(17)/pari(23+lst(32))
259  ENDIF
260 
261  ENDIF
262  RETURN
263  END