ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
tog4.F
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file tog4.F
1 *
2 * ********************************************************************
3 * * License and Disclaimer *
4 * * *
5 * * The Geant4 software is copyright of the Copyright Holders of *
6 * * the Geant4 Collaboration. It is provided under the terms and *
7 * * conditions of the Geant4 Software License, included in the file *
8 * * LICENSE and available at http://cern.ch/geant4/license . These *
9 * * include a list of copyright holders. *
10 * * *
11 * * Neither the authors of this software system, nor their employing *
12 * * institutes,nor the agencies providing financial support for this *
13 * * work make any representation or warranty, express or implied, *
14 * * regarding this software system or assume any liability for its *
15 * * use. Please see the license in the file LICENSE and URL above *
16 * * for the full disclaimer and the limitation of liability. *
17 * * *
18 * * This code implementation is the result of the scientific and *
19 * * technical work of the GEANT4 collaboration. *
20 * * By using, copying, modifying or distributing the software (or *
21 * * any work based on the software) you agree to acknowledge its *
22 * * use in resulting scientific publications, and indicate your *
23 * * acceptance of all terms of the Geant4 Software license. *
24 * ********************************************************************
25 *
26 *
27 *
28  subroutine tog4
29 ************************************************************************
30 *
31 * tog4
32 *
33 * Perform the translation to G4
34 *
35 ************************************************************************
36  implicit none
37 #include "gcbank.inc"
38  integer maxdivols
39  parameter(maxdivols=20000)
40  integer nvol, nrotm, nmate, ntmed, nset, i, jma, nmixt, k, nin,
41  > jdiv, jd, iaxis, ivo, ndiv, numed, npar, natt, ivol, jin,
42  > nparv, npard, nr, irot, konly, nwbuf, isvol, nmat, ifield,
43  > nbits(5000), idtyp, nwhi, nwdi, iset, idet, j, in, jmx,
44  > jdh, jdd, jdu, ndet, nn, nupar, npos, ndvol, ndivols, ii,
45  > npositioned, iia(10000), imate, smixt
46 
47  real c0, step, x, y, z, a, dens, radl, absl, fact(5000),
48  > fieldm, tmaxfd, stemax, deemax, epsil, stmin, orig(5000),
49  > upar(5000)
50 
51  character shape*4, name*4, dname*4, chonly*4, chmat*20, chtmed*20,
52  > chset*4, chdet*4, chnms(5000)*4, divols(maxdivols)*4
53 *
54  npositioned = 0
55 *
56 *** count materials and convert
57  call bankcnt(jmate,iia, nmate)
58  print *,'Materials: ',nmate
59  do imate=1,nmate
60  ii=iia(imate)
61  jma = lq(jmate-ii)
62  call uhtoc(iq(jma+1),4,chmat,20)
63  a = q(jma+6)
64  z = q(jma+7)
65  dens = q(jma+8)
66  radl = q(jma+9)
67  absl = q(jma+10)
68  nwbuf = iq(jma-1)-11
69  if (jma.gt.0) then
70  smixt=q(jma+11)
71  nmixt=abs(smixt)
72  if (nmixt.le.1) then
73  write(6,101) imate, chmat, a, z, dens, radl, absl
74  call ksmate(ii, chmat, a, z, dens, radl, absl,
75  > q(jma+12), nwbuf)
76  else
77  jmx = lq(jma-5)
78  write(6,102) imate, chmat, a, z, dens, radl, absl,
79  > (j,q(jmx+j),q(jmx+nmixt+j),q(jmx+2*nmixt+j),
80  > j=1,nmixt)
81  call ksmixt(ii, chmat, q(jmx+1), q(jmx+nmixt+1),
82  > dens, smixt, q(jmx+2*nmixt+1))
83  end if
84  end if
85  enddo
86  101 format(1x,i5,1x,a12,f6.2,f5.1,f8.2,2f9.2)
87  102 format(1x,i5,1x,a12,f6.2,f5.1,f8.2,2f9.2,1x,i2, f6.2, f5.1,
88  > f6.2/(57x, i2, f6.2, f5.1, f6.2))
89 *
90 *** count tracking media and convert
91  call bankcnt(jtmed,iia, ntmed)
92  print *,'Media: ',ntmed
93  do i=1,ntmed
94  ii=iia(i)
95  j = lq(jtmed-ii)
96  call uhtoc(iq(j+1),4,chtmed,20)
97  nmat = q(j+6)
98  isvol = q(j+7)
99  ifield = q(j+8)
100  fieldm = q(j+9)
101  tmaxfd = q(j+10)
102  stemax = q(j+11)
103  deemax = q(j+12)
104  epsil = q(j+13)
105  stmin = q(j+14)
106  nwbuf = iq(j-1) -14
107  call kstmed(ii,chtmed,nmat,isvol,ifield,fieldm,tmaxfd,stemax,
108  + deemax,epsil,stmin,q(j+15),nwbuf)
109  enddo
110 *
111 *** count rotation matrices and convert
112  call bankcnt(jrotm,iia, nrotm)
113  print *,'Rotms: ',nrotm
114  do i=1,nrotm
115  ii=iia(i)
116  j = lq(jrotm-ii)
117  call ksrotm(ii,q(j+11),q(j+12),q(j+13),q(j+14),q(j+15),q(j+16))
118  enddo
119 *
120 *** count volumes
121  npos = 0
122  call bankcnt(jvolum,iia, nvol)
123  print *,'Volumes: ',nvol
124 *** pull out the names of the volumes which are subvolumes of
125 *** divided volumes (gsvolu should not be called on these)
126  ndivols = 0
127  do i=1, nvol
128  ii=iia(i)
129  j = lq(jvolum-ii)
130  nin = q(j+3)
131  if (nin.lt.0) then
132  jdiv = lq(j-1)
133  ivo = q(jdiv+2)
134  call uhtoc(iq(jvolum+ivo),4,dname,4)
135  ndivols = ndivols +1
136  if (ndivols.gt.maxdivols) then
137  ndivols = maxdivols
138  print *,
139  + '!!!ERROR!!! ndivols array exhausted. ',
140  + 'Too many divisions.'
141  endif
142  divols(ndivols) = dname
143  endif
144  enddo
145 *** create the logical volumes (gsvolu)
146  ndvol = 0
147  do i=1, nvol
148  ii=iia(i)
149  j = lq(jvolum-ii)
150  call uhtoc(iq(jvolum+ii),4,name,4)
151  call jshape(q(j+2),shape)
152  nin = q(j+3)
153  numed = q(j+4)
154  npar = q(j+5)
155  natt = q(j+6)
156  do k=1, ndivols
157  if (divols(k).eq.name) then
158  ndvol = ndvol +1
159 c print *,'Division volume ',name,'; no gsvolu call.'
160  goto 11
161  endif
162  enddo
163  call ksvolu(name, shape, numed, q(j+7), npar, ivol)
164  11 continue
165  enddo
166  print *,'Divided volumes: ',ndvol
167 
168 *** properties of the mother volume
169  call uhtoc(iq(jvolum+1),4,name,4)
170  j=lq(jvolum-1)
171  call jshape(q(j+2),shape)
172  print *,'mother volume: ',name,' shape: ',shape
173 *** convert physical volumes
174  do i=1, nvol
175  ii=iia(i)
176  j = lq(jvolum-ii)
177  call uhtoc(iq(jvolum+ii),4,name,4)
178  nin = q(j+3)
179  numed = q(j+4)
180  npar = q(j+5)
181  if (nin.lt.0) then
182 * ! divided volume
183  jdiv = lq(j-1)
184  iaxis = q(jdiv+1)
185  ivo = q(jdiv+2)
186  call uhtoc(iq(jvolum+ivo),4,dname,4)
187  jd = lq(jvolum-ivo)
188  numed = q(jd+4)
189  ndiv = q(jdiv+3)
190  c0 = q(jdiv+4)
191  step = q(jdiv+5)
192  call ksdvn2(dname, name, ndiv, iaxis, c0, numed)
193  else if (nin.gt.0) then
194 * ! volume not divided. Handle positioning of daughter vols
195  do in=1, nin
196  jin = lq(j-in)
197  ivo = q(jin+2)
198  call uhtoc(iq(jvolum+ivo),4,dname,4)
199  jd = lq(jvolum-ivo)
200  nparv = q(jd+5) ! NPAR declared in the GSVOLU call
201  nr = q(jin+3)
202  irot = q(jin+4)
203  x = q(jin+5)
204  y = q(jin+6)
205  z = q(jin+7)
206  konly = q(jin+8)
207  if (konly.ne.0) then
208  chonly = 'ONLY'
209  else
210  chonly = 'MANY'
211  endif
212  npard = iq(jin-1) -9
213  npositioned = npositioned +1
214  if (nparv.eq.0) then
215 * ! use GSPOSP
216  call ksposp(dname, nr, name, x, y, z, irot, chonly,
217  + q(jin+10), npard)
218  else
219 * ! GSPOS
220  call kspos(dname, nr, name, x, y, z, irot, chonly)
221  endif
222  enddo
223  endif
224  enddo
225 *
226 *** count sensitive detectors
227  call bankcnt(jset,iia, nset)
228  print *,'Sets: ',nset
229  do i=1,nset
230  ii=iia(i)
231  j = lq(jset-ii)
232  call uhtoc(iq(jset+i),4,chset,4)
233  ndet = iq(j-1)
234  do k=1,ndet
235  jd = lq(j-k)
236  call uhtoc(iq(j+k),4,chdet,4)
237  call gfdet(chset, chdet, nn, chnms, nbits, idtyp,
238  + nwhi, nwdi, iset, idet)
239  call ksdet(chset, chdet, nn, chnms, nbits, idtyp,
240  + nwhi, nwdi, iset, idet)
241  jdh = lq(jd-1)
242  if (jdh.ne.0) then
243  call gfdeth(chset,chdet,nn,chnms,nbits,orig,fact)
244  call ksdeth(chset,chdet,nn,chnms,nbits,orig,fact)
245  endif
246  jdd = lq(jd-2)
247  if (jdd.ne.0) then
248  call gfdetd(chset,chdet,nn,chnms,nbits)
249  call ksdetd(chset,chdet,nn,chnms,nbits)
250  endif
251  jdu = lq(jd-3)
252  if (jdu.ne.0) then
253  call gfdetu(chset,chdet,100,nupar,upar)
254  call ksdetu(chset,chdet,nupar,upar)
255  endif
256  enddo
257  enddo
258  print *,'Positioned volumes (gspos, gsposp):',npositioned
259 *
260  call kgclos
261 *
262  end
263 
264  subroutine bankcnt(link,iia,nbanks)
265 ************************************************************************
266 ************************************************************************
267  implicit none
268 #include "gcbank.inc"
269  integer i, link, nbanks, iia(*)
270 *
271  nbanks=0
272  if (link.eq.0) return
273 C* do i=1,9999999
274  do i=1,iq(link-2)
275 C* if(lq(link-nbanks-1).eq.0.or.iq(link-2).eq.nbanks) goto 10
276  if(lq(link-i).ne.0)then
277  nbanks = nbanks +1
278  iia(nbanks)=i
279  endif
280  enddo
281  10 continue
282  end