Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
pyhikmap.f
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file pyhikmap.f
1 
2 C*********************************************************************
3 
4  SUBROUTINE pyhikmap(IVAR,MVAR,VVAR)
5 
6 C...Maps a uniform distribution into a distribution of a kinematical
7 C...variable according to one of the possibilities allowed. It is
8 C...assumed that kinematical limits have been set by a PYKLIM call.
9  common/ludat2/kchg(500,3),pmas(500,4),parf(2000),vckm(4,4)
10  SAVE /ludat2/
11  common/pyhiint1/mint(400),vint(400)
12  SAVE /pyhiint1/
13  common/pyhiint2/iset(200),kfpr(200,2),coef(200,20),icol(40,4,2)
14  SAVE /pyhiint2/
15 
16 C...Convert VVAR to tau variable.
17  isub=mint(1)
18  IF(ivar.EQ.1) THEN
19  taumin=vint(11)
20  taumax=vint(31)
21  IF(mvar.EQ.3.OR.mvar.EQ.4) THEN
22  taure=vint(73)
23  gamre=vint(74)
24  ELSEIF(mvar.EQ.5.OR.mvar.EQ.6) THEN
25  taure=vint(75)
26  gamre=vint(76)
27  ENDIF
28  IF(mint(43).EQ.1.AND.(iset(isub).EQ.1.OR.iset(isub).EQ.2)) THEN
29  tau=1.
30  ELSEIF(mvar.EQ.1) THEN
31  tau=taumin*(taumax/taumin)**vvar
32  ELSEIF(mvar.EQ.2) THEN
33  tau=taumax*taumin/(taumin+(taumax-taumin)*vvar)
34  ELSEIF(mvar.EQ.3.OR.mvar.EQ.5) THEN
35  ratgen=(taure+taumax)/(taure+taumin)*taumin/taumax
36  tau=taure*taumin/((taure+taumin)*ratgen**vvar-taumin)
37  ELSE
38  aupp=atan((taumax-taure)/gamre)
39  alow=atan((taumin-taure)/gamre)
40  tau=taure+gamre*tan(alow+(aupp-alow)*vvar)
41  ENDIF
42  vint(21)=min(taumax,max(taumin,tau))
43 
44 C...Convert VVAR to y* variable.
45  ELSEIF(ivar.EQ.2) THEN
46  ystmin=vint(12)
47  ystmax=vint(32)
48  IF(mint(43).EQ.1) THEN
49  yst=0.
50  ELSEIF(mint(43).EQ.2) THEN
51  IF(iset(isub).LE.2) yst=-0.5*log(vint(21))
52  IF(iset(isub).GE.3) yst=-0.5*log(vint(26))
53  ELSEIF(mint(43).EQ.3) THEN
54  IF(iset(isub).LE.2) yst=0.5*log(vint(21))
55  IF(iset(isub).GE.3) yst=0.5*log(vint(26))
56  ELSEIF(mvar.EQ.1) THEN
57  yst=ystmin+(ystmax-ystmin)*sqrt(vvar)
58  ELSEIF(mvar.EQ.2) THEN
59  yst=ystmax-(ystmax-ystmin)*sqrt(1.-vvar)
60  ELSE
61  aupp=atan(exp(ystmax))
62  alow=atan(exp(ystmin))
63  yst=log(tan(alow+(aupp-alow)*vvar))
64  ENDIF
65  vint(22)=min(ystmax,max(ystmin,yst))
66 
67 C...Convert VVAR to cos(theta-hat) variable.
68  ELSEIF(ivar.EQ.3) THEN
69  rm34=2.*vint(63)*vint(64)/(vint(21)*vint(2))**2
70  rsqm=1.+rm34
71  IF(2.*vint(71)**2/(vint(21)*vint(2)).LT.0.0001) rm34=max(rm34,
72  & 2.*vint(71)**2/(vint(21)*vint(2)))
73  ctnmin=vint(13)
74  ctnmax=vint(33)
75  ctpmin=vint(14)
76  ctpmax=vint(34)
77  IF(mvar.EQ.1) THEN
78  aneg=ctnmax-ctnmin
79  apos=ctpmax-ctpmin
80  IF(aneg.GT.0..AND.vvar*(aneg+apos).LE.aneg) THEN
81  vctn=vvar*(aneg+apos)/aneg
82  cth=ctnmin+(ctnmax-ctnmin)*vctn
83  ELSE
84  vctp=(vvar*(aneg+apos)-aneg)/apos
85  cth=ctpmin+(ctpmax-ctpmin)*vctp
86  ENDIF
87  ELSEIF(mvar.EQ.2) THEN
88  rmnmin=max(rm34,rsqm-ctnmin)
89  rmnmax=max(rm34,rsqm-ctnmax)
90  rmpmin=max(rm34,rsqm-ctpmin)
91  rmpmax=max(rm34,rsqm-ctpmax)
92  aneg=log(rmnmin/rmnmax)
93  apos=log(rmpmin/rmpmax)
94  IF(aneg.GT.0..AND.vvar*(aneg+apos).LE.aneg) THEN
95  vctn=vvar*(aneg+apos)/aneg
96  cth=rsqm-rmnmin*(rmnmax/rmnmin)**vctn
97  ELSE
98  vctp=(vvar*(aneg+apos)-aneg)/apos
99  cth=rsqm-rmpmin*(rmpmax/rmpmin)**vctp
100  ENDIF
101  ELSEIF(mvar.EQ.3) THEN
102  rmnmin=max(rm34,rsqm+ctnmin)
103  rmnmax=max(rm34,rsqm+ctnmax)
104  rmpmin=max(rm34,rsqm+ctpmin)
105  rmpmax=max(rm34,rsqm+ctpmax)
106  aneg=log(rmnmax/rmnmin)
107  apos=log(rmpmax/rmpmin)
108  IF(aneg.GT.0..AND.vvar*(aneg+apos).LE.aneg) THEN
109  vctn=vvar*(aneg+apos)/aneg
110  cth=rmnmin*(rmnmax/rmnmin)**vctn-rsqm
111  ELSE
112  vctp=(vvar*(aneg+apos)-aneg)/apos
113  cth=rmpmin*(rmpmax/rmpmin)**vctp-rsqm
114  ENDIF
115  ELSEIF(mvar.EQ.4) THEN
116  rmnmin=max(rm34,rsqm-ctnmin)
117  rmnmax=max(rm34,rsqm-ctnmax)
118  rmpmin=max(rm34,rsqm-ctpmin)
119  rmpmax=max(rm34,rsqm-ctpmax)
120  aneg=1./rmnmax-1./rmnmin
121  apos=1./rmpmax-1./rmpmin
122  IF(aneg.GT.0..AND.vvar*(aneg+apos).LE.aneg) THEN
123  vctn=vvar*(aneg+apos)/aneg
124  cth=rsqm-1./(1./rmnmin+aneg*vctn)
125  ELSE
126  vctp=(vvar*(aneg+apos)-aneg)/apos
127  cth=rsqm-1./(1./rmpmin+apos*vctp)
128  ENDIF
129  ELSEIF(mvar.EQ.5) THEN
130  rmnmin=max(rm34,rsqm+ctnmin)
131  rmnmax=max(rm34,rsqm+ctnmax)
132  rmpmin=max(rm34,rsqm+ctpmin)
133  rmpmax=max(rm34,rsqm+ctpmax)
134  aneg=1./rmnmin-1./rmnmax
135  apos=1./rmpmin-1./rmpmax
136  IF(aneg.GT.0..AND.vvar*(aneg+apos).LE.aneg) THEN
137  vctn=vvar*(aneg+apos)/aneg
138  cth=1./(1./rmnmin-aneg*vctn)-rsqm
139  ELSE
140  vctp=(vvar*(aneg+apos)-aneg)/apos
141  cth=1./(1./rmpmin-apos*vctp)-rsqm
142  ENDIF
143  ENDIF
144  IF(cth.LT.0.) cth=min(ctnmax,max(ctnmin,cth))
145  IF(cth.GT.0.) cth=min(ctpmax,max(ctpmin,cth))
146  vint(23)=cth
147 
148 C...Convert VVAR to tau' variable.
149  ELSEIF(ivar.EQ.4) THEN
150  tau=vint(11)
151  taupmn=vint(16)
152  taupmx=vint(36)
153  IF(mint(43).EQ.1) THEN
154  taup=1.
155  ELSEIF(mvar.EQ.1) THEN
156  taup=taupmn*(taupmx/taupmn)**vvar
157  ELSE
158  aupp=(1.-tau/taupmx)**4
159  alow=(1.-tau/taupmn)**4
160  taup=tau/(1.-(alow+(aupp-alow)*vvar)**0.25)
161  ENDIF
162  vint(26)=min(taupmx,max(taupmn,taup))
163  ENDIF
164 
165  RETURN
166  END