ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
g3routines.F
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file g3routines.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 #define CALL_GEANT
29 
30 #ifndef CALL_GEANT
31  subroutine gsvolu(name, shape, nmed, par, npar, ivol)
32 #else
33  subroutine ksvolu(name, shape, nmed, par, npar, ivol)
34 #endif
35 ************************************************************************
36 ************************************************************************
37  implicit none
38  character name*4, shape*4, fmt*150
39  integer nmed, npar, ivol, k
40  real par(npar)
41  character rname*6
42 #include "G3toG4.inc"
43  data rname /'GSVOLU'/
44 *
45  call check_lines
46 #ifdef CALL_GEANT
47  if (dogeom) call gsvolu(name, shape, nmed, par, npar, ivol)
48 #endif
49  if (npar.ne.0) call checkshape(name, shape, par, npar)
50 *
51  if (lunlist.ne.0) then
52 * write(lunlist,
53 * + '(a4,1x,a6,1x,a4,1x,a4,2i5,<npar>e15.8)')
54 * + context, rname, name, shape, nmed, npar,
55 * + (par(k),k=1,npar)
56  write(fmt,'(A,I2,A)')'(a4,1x,a6,1x,a4,1x,a4,2i5,',max(npar,1),
57  > '(1x,e16.8))'
58  write(lunlist,fmt) context, rname, name, shape, nmed, npar,
59  + (par(k),k=1,npar)
60  endif
61  if (luncode.ne.0) then
62  write(luncode,'(''{'')')
63  call g3ldpar(par,npar)
64  write(luncode,1000) name, shape, nmed, npar
65  1000 format('G4gsvolu(name="',a,'",shape="',a,'",nmed=',i5,
66  + ',par,npar=',i4,');')
67  write(luncode,'(''}'')')
68  endif
69 *
70  end
71 *
72 #ifndef CALL_GEANT
73  subroutine gspos(name, num, moth, x, y, z, irot, only)
74 #else
75  subroutine kspos(name, num, moth, x, y, z, irot, only)
76 #endif
77 ************************************************************************
78 ************************************************************************
79  implicit none
80  character name*4, moth*4, only*4
81  integer num, irot
82  real x, y, z
83  character rname*6
84 #include "G3toG4.inc"
85  data rname /'GSPOS '/
86 *
87  call check_lines
88 #ifdef CALL_GEANT
89  if (dogeom) call gspos(name, num, moth, x, y, z, irot, only)
90 #endif
91  if (lunlist.ne.0) then
92  write(lunlist,
93  + '(a4,1x,a6,1x,a4,i5,1x,a4,3(1x,e16.8),i5,1x,a4)')
94  + context, rname, name, num, moth, x, y, z, irot, only
95  endif
96  if (luncode.ne.0) then
97  write(luncode,'(''{'')')
98  call rtocp('x',x)
99  call rtocp('y',y)
100  call rtocp('z',z)
101  write(luncode,1000) name,num,moth,irot,only
102  1000 format('G4gspos(name="',a,'",num=',i5,',moth="',a,
103  + '",x,y,z,irot=',i5,',only="',a,'");')
104  write(luncode,'(''}'')')
105  endif
106 *
107  end
108 *
109 #ifndef CALL_GEANT
110  subroutine gsposp(name, num, moth, x, y, z, irot, only, par, npar)
111 #else
112  subroutine ksposp(name, num, moth, x, y, z, irot, only, par, npar)
113 #endif
114 ************************************************************************
115 ************************************************************************
116  implicit none
117  character name*4, moth*4, only*4
118  integer num, irot, npar, k
119  real x, y, z, par(npar)
120  character rname*6, fmt*150
121 #include "G3toG4.inc"
122  data rname /'GSPOSP'/
123 *
124  call check_lines
125 #ifdef CALL_GEANT
126  if (dogeom) call gsposp(name, num, moth, x, y, z, irot, only,
127  + par, npar)
128 #endif
129  if (lunlist.ne.0) then
130  do k=1,npar
131  if (abs(par(k)).gt.1.e10) then
132  print *,'Warning: huge junk value in PAR for GSPOS'
133  print *,' zeroed out. Volume ',name
134  par(k) = 0.
135  endif
136  enddo
137 * write(lunlist,
138 * + '(a4,1x,a6,1x,a4,i5,1x,a4,3e15.8,i5,1x,a4,
139 * + i5,<npar>e15.8)')
140 * + context, rname, name, num, moth, x, y, z, irot, only,
141 * + npar,
142 * + (par(k),k=1,npar)
143  write(fmt,'(A,A,I2,A)')
144  > '(a4,1x,a6,1x,a4,i5,1x,a4,3(1x,e16.8),',
145  + 'i5,1x,a4,i5,',max(npar,1),'(1x,e16.8))'
146  write(lunlist,fmt)
147  + context, rname, name, num, moth, x, y, z, irot, only,
148  + npar,
149  + (par(k),k=1,npar)
150  endif
151  if (luncode.ne.0) then
152  write(luncode,'(''{'')')
153  call rtocp('x',x)
154  call rtocp('y',y)
155  call rtocp('z',z)
156  call g3ldpar(par,npar)
157  write(luncode,1000) name,num,moth,irot,only,npar
158  1000 format('G4gsposp(name="',a,'",num=',i5,',moth="',a,
159  + '",x,y,z,irot=',i5,',only="',a,'",par,npar=',i4,');')
160  write(luncode,'(''}'')')
161  endif
162 *
163  end
164 *
165 #ifndef CALL_GEANT
166  subroutine gsatt(name, attr, ival)
167 #else
168  subroutine ksatt(name, attr, ival)
169 #endif
170 ************************************************************************
171 ************************************************************************
172  implicit none
173  character name*4, attr*4
174  integer ival
175  character rname*6
176 #include "G3toG4.inc"
177  data rname /'GSATT '/
178 *
179  call check_lines
180 #ifdef CALL_GEANT
181  if (dogeom) call gsatt(name, attr, ival)
182 #endif
183  if (lunlist.ne.0) then
184  write(lunlist,
185  + '(a4,1x,a6,1x,a4,1x,a4,i12)')
186  + context, rname, name, attr, ival
187  endif
188  if (luncode.ne.0) then
189  write(luncode,'(''{'')')
190  write(luncode,1000) name,attr,ival
191  1000 format('G4gsatt(name="',a,'",attr="',a,'",ival=',i10,');')
192  write(luncode,'(''}'')')
193  endif
194 *
195  end
196 *
197 #ifndef CALL_GEANT
198  subroutine gsrotm(irot, theta1, phi1, theta2, phi2,
199  + theta3, phi3)
200 #else
201  subroutine ksrotm(irot, theta1, phi1, theta2, phi2,
202  + theta3, phi3)
203 #endif
204 ************************************************************************
205 ************************************************************************
206  implicit none
207  integer irot
208  real theta1, phi1, theta2, phi2, theta3, phi3
209  character rname*6
210 #include "G3toG4.inc"
211  data rname /'GSROTM'/
212 *
213  call check_lines
214 #ifdef CALL_GEANT
215  if (dogeom) call gsrotm(irot, theta1, phi1, theta2, phi2,
216  + theta3, phi3)
217 #endif
218  if (lunlist.ne.0) then
219  write(lunlist,
220  + '(a4,1x,a6,i5,6f11.5)')
221  + context, rname, irot, theta1, phi1, theta2, phi2,
222  + theta3, phi3
223  endif
224  if (luncode.ne.0) then
225  write(luncode,'(''{'')')
226  call rtocp('theta1',theta1)
227  call rtocp('phi1',phi1)
228  call rtocp('theta2',theta2)
229  call rtocp('phi2',phi2)
230  call rtocp('theta3',theta3)
231  call rtocp('phi3',phi3)
232  write(luncode,1000) irot
233  1000 format('G4gsrotm(irot=',i5,
234  + ',theta1,phi1,theta2,phi2,theta3,phi3);')
235  write(luncode,'(''}'')')
236  endif
237 *
238  end
239 *
240 #ifndef CALL_GEANT
241  subroutine gsdvn(name, moth, ndiv, iaxis)
242 #else
243  subroutine ksdvn(name, moth, ndiv, iaxis)
244 #endif
245 ************************************************************************
246 ************************************************************************
247  implicit none
248  character name*4, moth*4
249  integer ndiv, iaxis
250  character rname*6
251 #include "G3toG4.inc"
252  data rname /'GSDVN '/
253 *
254  call check_lines
255 #ifdef CALL_GEANT
256  if (dogeom) call gsdvn(name, moth, ndiv, iaxis)
257 #endif
258  if (lunlist.ne.0) then
259  write(lunlist,
260  + '(a4,1x,a6,1x,a4,1x,a4,i5,i3)')
261  + context, rname, name, moth, ndiv, iaxis
262  endif
263  if (luncode.ne.0) then
264  write(luncode,'(''{'')')
265  write(luncode,1000) name, moth, ndiv, iaxis
266  1000 format('G4gsdvn(name="',a,'",moth="',a,'",ndiv=',i3,
267  + ',iaxis=',i1,');')
268  write(luncode,'(''}'')')
269  endif
270 *
271  end
272 *
273 #ifndef CALL_GEANT
274  subroutine gsdvt(name, moth, step, iaxis, numed, ndvmx)
275 #else
276  subroutine ksdvt(name, moth, step, iaxis, numed, ndvmx)
277 #endif
278 ************************************************************************
279 ************************************************************************
280  implicit none
281  character name*4, moth*4
282  real step
283  integer iaxis, numed, ndvmx
284  character rname*6
285 #include "G3toG4.inc"
286  data rname /'GSDVT '/
287 *
288  call check_lines
289 #ifdef CALL_GEANT
290  if (dogeom) call gsdvt(name, moth, step, iaxis, numed, ndvmx)
291 #endif
292  if (lunlist.ne.0) then
293  write(lunlist,
294  + '(a4,1x,a6,1x,a4,1x,a4,(1x,e16.8),3i5)')
295  + context, rname, name, moth, step, iaxis, numed, ndvmx
296  endif
297  if (luncode.ne.0) then
298  write(luncode,'(''{'')')
299  call rtocp('step',step)
300  write(luncode,1000) name,moth,iaxis,numed,ndvmx
301  1000 format('G4gsdvt(name="',a,'",moth="',a,'",step,iaxis=',
302  + i1,',numed=',i4,',ndvmx=',i4,');')
303  write(luncode,'(''}'')')
304  endif
305 *
306  end
307 *
308 #ifndef CALL_GEANT
309  subroutine gsdvx(name, moth, ndiv, iaxis, step, c0, numed, ndvmx)
310 #else
311  subroutine ksdvx(name, moth, ndiv, iaxis, step, c0, numed, ndvmx)
312 #endif
313 ************************************************************************
314 ************************************************************************
315  implicit none
316  character name*4, moth*4
317  integer ndiv, iaxis, numed, ndvmx
318  real step, c0
319  character rname*6
320 #include "G3toG4.inc"
321  data rname /'GSDVX '/
322 *
323  call check_lines
324 #ifdef CALL_GEANT
325  if (dogeom) call gsdvx(name, moth, ndiv, iaxis, step, c0, numed,
326  + ndvmx)
327 #endif
328  if (lunlist.ne.0) then
329  write(lunlist,
330  + '(a4,1x,a6,1x,a4,1x,a4,i5,i3,2(1x,e16.8),2i5)')
331  + context, rname, name, moth, ndiv, iaxis,step, c0,
332  + numed, ndvmx
333  endif
334  if (luncode.ne.0) then
335  write(luncode,'(''{'')')
336  call rtocp('step',step)
337  call rtocp('c0',c0)
338  write(luncode,1000) name,moth,ndiv,iaxis,numed,ndvmx
339  1000 format('G4gsdvx(name="',a,'",moth="',a,'",ndiv=',i3,',iaxis=',
340  + i1,',step,c0,numed=',i4,',ndvmx=',i4,');')
341  write(luncode,'(''}'')')
342  endif
343 *
344  end
345 *
346 #ifndef CALL_GEANT
347  subroutine gsdvn2(name, moth, ndiv, iaxis, c0, numed)
348 #else
349  subroutine ksdvn2(name, moth, ndiv, iaxis, c0, numed)
350 #endif
351 ************************************************************************
352 ************************************************************************
353  implicit none
354  character name*4, moth*4
355  integer ndiv, iaxis, numed
356  real c0
357  character rname*6
358 #include "G3toG4.inc"
359  data rname /'GSDVN2'/
360 *
361  call check_lines
362 #ifdef CALL_GEANT
363  if (dogeom) call gsdvn2(name, moth, ndiv, iaxis, c0, numed)
364 #endif
365  if (lunlist.ne.0) then
366  write(lunlist,
367  + '(a4,1x,a6,1x,a4,1x,a4,i5,i3,(1x,e16.8),i5)')
368  + context, rname, name, moth, ndiv, iaxis, c0, numed
369  endif
370  if (luncode.ne.0) then
371  write(luncode,'(''{'')')
372  call rtocp('c0',c0)
373  write(luncode, 1000) name,moth,ndiv,iaxis,numed
374  1000 format('G4gsdvn2(name="',a,'",moth="',a,'",ndiv=',i3,',iaxis=',
375  + i1,',c0,numed=',i4,');')
376  write(luncode,'(''}'')')
377  endif
378 *
379  end
380 *
381 #ifndef CALL_GEANT
382  subroutine gsdvt2(name, moth, step, iaxis, c0, numed, ndvmx)
383 #else
384  subroutine ksdvt2(name, moth, step, iaxis, c0, numed, ndvmx)
385 #endif
386 ************************************************************************
387 ************************************************************************
388  implicit none
389  character name*4, moth*4
390  integer iaxis, numed, ndvmx
391  real step, c0
392  character rname*6
393 #include "G3toG4.inc"
394  data rname /'GSDVT2'/
395 *
396  call check_lines
397 #ifdef CALL_GEANT
398  if (dogeom) call gsdvt2(name, moth, step, iaxis, c0, numed, ndvmx)
399 #endif
400  if (lunlist.ne.0) then
401  write(lunlist,
402  + '(a4,1x,a6,1x,a4,1x,a4,(1x,e16.8),i3,(1x,e16.8),2i5)')
403  + context, rname, name, moth, step, iaxis, c0, numed, ndvmx
404  endif
405  if (luncode.ne.0) then
406  write(luncode,'(''{'')')
407  call rtocp('step',step)
408  call rtocp('c0',c0)
409  write(luncode,1000) name,moth,iaxis,numed,ndvmx
410  1000 format('G4gsdvt2(name="',a,'",moth="',a,'",step,iaxis=',
411  + i1,',c0,numed=',i4,',ndvmx=',i4,');')
412  write(luncode,'(''}'')')
413  endif
414 *
415  end
416 *
417 #ifndef CALL_GEANT
418  subroutine gsmate(imate, name, a, z, dens, radl, absl, ubf, nwbf)
419 #else
420  subroutine ksmate(imate, name, a, z, dens, radl, absl, ubf, nwbf)
421 #endif
422 ************************************************************************
423 ************************************************************************
424  implicit none
425  character name*(*)
426  integer imate, nwbf, k
427  real a, z, dens, radl, absl, ubf(nwbf)
428  character rname*6, fmt*150
429 #include "G3toG4.inc"
430  data rname /'GSMATE'/
431 *
432  call check_lines
433 #ifdef CALL_GEANT
434  if (dogeom) call gsmate
435  + (imate, name, a, z, dens, radl, absl, ubf, nwbf)
436 #endif
437  if (lunlist.ne.0) then
438  write(fmt,'(A,I3,A)')
439  > '(a4,1x,a6,i5,1x,''"'',a,''"'',4(1x,e16.8),i3,',
440  > max(nwbf,1),'(1x,e16.8))'
441  write(lunlist,fmt)
442  + context, rname, imate, name, a, z, dens, radl,
443  + nwbf, (ubf(k), k=1,nwbf)
444  endif
445  if (luncode.ne.0) then
446  write(luncode,'(''{'')')
447  call rtocp('a',a)
448  call rtocp('z',z)
449  call rtocp('dens',dens)
450  call rtocp('radl',radl)
451  call g3ldpar(ubf,nwbf)
452  write(luncode,1000) imate, name, nwbf
453  1000 format('G4gsmate(imate=',i4,',name="',a,
454  + '",a,z,dens,radl,npar=',i4,',par);')
455  write(luncode,'(''}'')')
456  endif
457 *
458  end
459 *
460 #ifndef CALL_GEANT
461  subroutine gsmixt(imate, name, a, z, dens, nlmat, wmat)
462 #else
463  subroutine ksmixt(imate, name, a, z, dens, nlmat, wmat)
464 #endif
465 ************************************************************************
466 ************************************************************************
467  implicit none
468  character name*(*)
469  integer imate, nlmat, k, nlmata
470  real a(*), z(*), dens, wmat(*)
471  character rname*6, fmt*150
472 #include "G3toG4.inc"
473  data rname /'GSMIXT'/
474 *
475  call check_lines
476 #ifdef CALL_GEANT
477  if (dogeom) call gsmixt
478  + (imate, name, a, z, dens, nlmat, wmat)
479 #endif
480  if (lunlist.ne.0) then
481  nlmata = abs(nlmat)
482  write(fmt,'(A,I3,A,I3,A,I3,A)')
483  + '(a4,1x,a6,i5,1x,''"'',a,''"'',1x,e16.8,1x,i3,',
484  > max(nlmata,1),
485  > '(1x,e16.8),',max(nlmata,1),'(1x,e16.8),',
486  > max(nlmata,1),'(1x,e16.8))'
487  write(lunlist,fmt)
488  + context, rname, imate, name, dens,
489  + nlmat,
490  + (a(k), k=1,abs(nlmat)),
491  + (z(k), k=1,abs(nlmat)),
492  + (wmat(k), k=1,abs(nlmat))
493  endif
494  if (luncode.ne.0) then
495  write(luncode,'(''{'')')
496  call rtocp('dens',dens)
497  call artocp('aa',a,abs(nlmat))
498  call artocp('zz',z,abs(nlmat))
499  call artocp('wmat',wmat,abs(nlmat))
500  write(luncode,1000) imate,name,nlmat
501  1000 format('G4gsmixt(imate=',i5,',name="',a,
502  + '",aa,zz,dens,nlmat=',i3,',wmat);')
503  write(luncode,'(''}'')')
504  endif
505 *
506  end
507 *
508 #ifndef CALL_GEANT
509  subroutine gstmed(
510  + itmed, name, nmat, isvol, ifield, fieldm,
511  + tmaxfd, stemax, deemax, epsil, stmin, ubuf, nwbuf)
512 #else
513  subroutine kstmed(
514  + itmed, name, nmat, isvol, ifield, fieldm,
515  + tmaxfd, stemax, deemax, epsil, stmin, ubuf, nwbuf)
516 #endif
517 ************************************************************************
518 ************************************************************************
519  implicit none
520  character name*(*)
521  integer itmed, nmat, isvol, ifield, nwbuf, k
522  real fieldm, tmaxfd, stemax, deemax, epsil, stmin, ubuf(nwbuf)
523  character rname*6, fmt*150
524 #include "G3toG4.inc"
525  data rname /'GSTMED'/
526 *
527  call check_lines
528 #ifdef CALL_GEANT
529  if (dogeom) call gstmed(
530  + itmed, name, nmat, isvol, ifield, fieldm,
531  + tmaxfd, stemax, deemax, epsil, stmin, ubuf, nwbuf)
532 #endif
533  if (lunlist.ne.0) then
534 * write(lunlist,
535 * + '(a4,1x,a6,i5,1x,''"'',a,''"'',3i3,6e15.8,i3,<nwbuf>e15.8)')
536 * + context, rname, itmed, name, nmat, isvol, ifield, fieldm,
537 * + tmaxfd, stemax, deemax, epsil, stmin,
538 * + nwbuf, (ubuf(k),k=1,nwbuf)
539  write(fmt,'(A,I3,A)')
540  > '(a4,1x,a6,i5,1x,''"'',a,''"'',3i3,6(1x,e16.8),i3,',
541  > max(nwbuf,1),'(1x,e16.8))'
542  write(lunlist,fmt)
543  + context, rname, itmed, name, nmat, isvol, ifield, fieldm,
544  + tmaxfd, stemax, deemax, epsil, stmin,
545  + nwbuf, (ubuf(k),k=1,nwbuf)
546  endif
547  if (luncode.ne.0) then
548  write(luncode,'(''{'')')
549  call rtocp('fieldm',fieldm)
550  call rtocp('tmaxfd',tmaxfd)
551  call rtocp('stemax',stemax)
552  call rtocp('deemax',deemax)
553  call rtocp('epsil',epsil)
554  call rtocp('stmin',stmin)
555  call g3ldpar(ubuf,nwbuf)
556  write(luncode,1000) itmed,name,nmat,isvol,ifield,nwbuf
557  1000 format('G4gstmed(itmed=',i4,',name="',a,'",nmat=',i4,
558  + ',isvol=',i2,',ifield=',i2,',',/
559  + ' fieldm,tmaxfd,stemax,deemax,epsil,stmin,par,npar=',
560  + i4,');')
561  write(luncode,'(''}'')')
562  endif
563 *
564  end
565 *
566 #ifndef CALL_GEANT
567  subroutine gstpar(itmed, chpar, parval)
568 #else
569  subroutine kstpar(itmed, chpar, parval)
570 #endif
571 ************************************************************************
572 ************************************************************************
573  implicit none
574  character chpar*(*)
575  integer itmed
576  real parval
577  character rname*6
578 #include "G3toG4.inc"
579  data rname /'GSTPAR'/
580 *
581  call check_lines
582 #ifdef CALL_GEANT
583  if (dogeom) call gstpar(itmed, chpar, parval)
584 #endif
585  if (lunlist.ne.0) then
586  write(lunlist,
587  + '(a4,1x,a6,i5,1x,a4,(1x,e16.8))')
588  + context, rname, itmed, chpar, parval
589  endif
590  if (luncode.ne.0) then
591  write(luncode,'(''{'')')
592  write(luncode,1000) itmed, chpar, parval
593  1000 format('G4gstpar(itmed=',i4,',chpar="',a,'",parval=',
594  + (1x,e16.8),');')
595  write(luncode,'(''}'')')
596  endif
597 *
598  end
599 *
600 #ifndef CALL_GEANT
601  subroutine gspart(
602  + ipart, chpar, itrtyp, amass, charge, tlife, ub, nwb)
603 #else
604  subroutine kspart(
605  + ipart, chpar, itrtyp, amass, charge, tlife, ub, nwb)
606 #endif
607 ************************************************************************
608 ************************************************************************
609  implicit none
610  character chpar*(*)
611  integer ipart, itrtyp, nwb, k
612  real amass, charge, tlife, ub(nwb)
613  character rname*6, fmt*150
614 #include "G3toG4.inc"
615  data rname /'GSPART'/
616 *
617  call check_lines
618 #ifdef CALL_GEANT
619  if (dogeom) call gspart(
620  + ipart, chpar, itrtyp, amass, charge, tlife, ub, nwb)
621 #endif
622  if (lunlist.ne.0) then
623 * write(lunlist,
624 * + '(a4,1x,a6,i5,1x,''"'',a,''"'',i3,3e15.8,i3,<nwb>e15.8)')
625 * + context, rname, ipart, chpar, itrtyp, amass, charge, tlife,
626 * + nwb, (ub(k), k=1,nwb)
627  write(fmt,'(A,I3,A)')
628  > '(a4,1x,a6,i5,1x,''"'',a,''"'',i3,3(1x,e16.8),i3,',
629  > max(nwb,1),'(1x,e16.8))'
630  write(lunlist,fmt)
631  + context, rname, ipart, chpar, itrtyp, amass, charge,
632  > tlife,
633  + nwb, (ub(k), k=1,nwb)
634  endif
635  if (luncode.ne.0) then
636  write(luncode,'(''{'')')
637  call rtocp('amass',amass)
638  call rtocp('charge',charge)
639  call rtocp('tlife',tlife)
640  call g3ldpar(ub,nwb)
641  write(luncode,1000) ipart,chpar,itrtyp,nwb
642  1000 format('G4gspart(ipart=',i8,',chpar="',a,'",itrtyp=',i8,
643  + ',amass,charge,'/' tlife,par,npar=',i4,');')
644  write(luncode,'(''}'')')
645  endif
646 *
647  end
648 *
649 #ifndef CALL_GEANT
650  subroutine gsdk(ipart, bratio, mode)
651 #else
652  subroutine ksdk(ipart, bratio, mode)
653 #endif
654 ************************************************************************
655 ************************************************************************
656  implicit none
657  integer ipart, mode(6)
658  real bratio(6)
659  character rname*6
660 #include "G3toG4.inc"
661  data rname /'GSDK '/
662 *
663  call check_lines
664 #ifdef CALL_GEANT
665  if (dogeom) call gsdk(ipart, bratio, mode)
666 #endif
667  if (lunlist.ne.0) then
668 *** 6 is prefixed to the arrays for consistency with other
669 *** array treatments (count precedes the array)
670  write(lunlist,
671  + '(a4,1x,a6,i5,i3,6(1x,e16.8),6i8)')
672  + context, rname, ipart, 6, bratio, mode
673  endif
674  if (luncode.ne.0) then
675  write(luncode,'(''{'')')
676  call artocp('bratio',bratio,6)
677  call aitocp('mode',mode,6)
678  write(luncode,1000) ipart
679  1000 format('G4gsdk(ipart=',i8,',bratio,mode);')
680  write(luncode,'(''}'')')
681  endif
682 *
683  end
684 *
685 #ifndef CALL_GEANT
686  subroutine gsdet(chset, chdet, nv, chnam, nbits, idtyp, nwhi,
687  + nwdi, iset, idet)
688 #else
689  subroutine ksdet(chset, chdet, nv, chnam, nbits, idtyp, nwhi,
690  + nwdi, iset, idet)
691 #endif
692 ************************************************************************
693 ************************************************************************
694  implicit none
695  integer nv, nbits(nv), idtyp, nwhi, nwdi, iset, idet, k
696  character rname*6, chset*4, chdet*4, chnam(nv)*4, fmt*150
697 #include "G3toG4.inc"
698  data rname /'GSDET '/
699 *
700  call check_lines
701 #ifdef CALL_GEANT
702  if (dogeom) call gsdet(chset, chdet, nv, chnam, nbits, idtyp,
703  + nwhi, nwdi, iset, idet)
704 #endif
705  if (lunlist.ne.0) then
706 * write(lunlist,
707 * + '(a4,1x,a6,1x,a4,1x,a4,i5,<nv>(1x,a4),<nv>i10,i10,2i5)')
708 * + context, rname, chset, chdet, nv, (chnam(k), k=1,nv),
709 * + (nbits(k), k=1,nv), idtyp, nwhi, nwdi
710  write(fmt,'(A,I3,A,I3,A)')'(a4,1x,a6,1x,a4,1x,a4,i5,',
711  > max(nv,1),'(1x,a4),',max(nv,1),'i10,i10,2i5)'
712  write(lunlist,fmt)
713  + context, rname, chset, chdet, nv, (chnam(k), k=1,nv),
714  + (nbits(k), k=1,nv), idtyp, nwhi, nwdi
715  endif
716  if (luncode.ne.0) then
717  write(luncode,'(''{'')')
718  call astocp('chnam',chnam,nv)
719  call aitocp('nbits',nbits,nv)
720  write(luncode,1000) chset, chdet, nv, idtyp, nwhi, nwdi
721  1000 format('G4gsdet(chset="',a,'",chdet="',a,'",nv=',i3,
722  + ',chnam,nbits,idtyp=',i8,','/
723  + ' nwhi=',i8,',nwdi=',i8,');')
724  write(luncode,'(''}'')')
725  endif
726 *
727  end
728 *
729 #ifndef CALL_GEANT
730  subroutine gsdetv(chset, chdet, idtyp, nwhi, nwdi, iset, idet)
731 #else
732  subroutine ksdetv(chset, chdet, idtyp, nwhi, nwdi, iset, idet)
733 #endif
734 ************************************************************************
735 ************************************************************************
736  implicit none
737  integer idtyp, nwhi, nwdi, iset, idet
738  character rname*6, chset*4, chdet*4
739 #include "G3toG4.inc"
740  data rname /'GSDETV'/
741 *
742  call check_lines
743 #ifdef CALL_GEANT
744  if (dogeom) call gsdetv(chset, chdet, idtyp,
745  + nwhi, nwdi, iset, idet)
746 #endif
747  if (lunlist.ne.0) then
748  write(lunlist,
749  + '(a4,1x,a6,1x,a4,1x,a4,i10,2i5)')
750  + context, rname, chset, chdet, idtyp, nwhi, nwdi
751  endif
752  if (luncode.ne.0) then
753  write(luncode,'(''{'')')
754  write(luncode,1000) chset, chdet, idtyp, nwhi, nwdi
755  1000 format('G4gsdetv(chset="',a,'",chdet="',a,'",idtyp=',i8,
756  + ',nwhi=',i8,',nwdi=',i8,');')
757  write(luncode,'(''}'')')
758  endif
759 *
760  end
761 *
762 #ifndef CALL_GEANT
763  subroutine gsdeta(chset, chdet, chali, nwhi, nwdi, iali)
764 #else
765  subroutine ksdeta(chset, chdet, chali, nwhi, nwdi, iali)
766 #endif
767 ************************************************************************
768 ************************************************************************
769  implicit none
770  integer nwhi, nwdi, iali
771  character rname*6, chset*4, chdet*4, chali*4
772 #include "G3toG4.inc"
773  data rname /'GSDETA'/
774 *
775  call check_lines
776 #ifdef CALL_GEANT
777  if (dogeom) call gsdeta(chset, chdet, chali, nwhi, nwdi, iali)
778 #endif
779  if (lunlist.ne.0) then
780  write(lunlist,
781  + '(a4,1x,a6,1x,a4,1x,a4,1x,a4,2i5)')
782  + context, rname, chset, chdet, chali, nwhi, nwdi
783  endif
784  if (luncode.ne.0) then
785  write(luncode,'(''{'')')
786  write(luncode,1000) chset, chdet, chali, nwhi, nwdi
787  1000 format('G4gsdeta(chset="',a,'",chdet="',a,'",chali="',a,
788  + '",nwhi=',i8,',nwdi=',i8,');')
789  write(luncode,'(''}'')')
790  endif
791 *
792  end
793 *
794 #ifndef CALL_GEANT
795  subroutine gsdeth(chset, chdet, nh, chnam, nbits, orig, fact)
796 #else
797  subroutine ksdeth(chset, chdet, nh, chnam, nbits, orig, fact)
798 #endif
799 ************************************************************************
800 ************************************************************************
801  implicit none
802  integer nh, nbits(nh), k
803  real orig(nh), fact(nh)
804  character rname*6, chset*4, chdet*4, chnam(nh)*4, fmt*150
805 #include "G3toG4.inc"
806  data rname /'GSDETH'/
807 *
808  call check_lines
809 #ifdef CALL_GEANT
810  if (dogeom) call gsdeth(chset, chdet, nh, chnam, nbits,
811  + orig, fact)
812 #endif
813  if (lunlist.ne.0) then
814 * write(lunlist,
815 * + '(a4,1x,a6,1x,a4,1x,a4,i5,<nh>(1x,a4),<nh>i5,<nh>e15.8,
816 * + <nh>e15.8)')
817 * + context, rname, chset, chdet, nh, (chnam(k), k=1,nh),
818 * + (nbits(k), k=1,nh), (orig(k), k=1,nh), (fact(k), k=1,nh)
819  write(fmt,'(A,I3,A,I3,A,I3,A,I3,A)')
820  > '(a4,1x,a6,1x,a4,1x,a4,i5,',max(nh,1),'(1x,a4),',
821  > max(nh,1),'i5,',max(nh,1),'(1x,e16.8),',max(nh,1),
822  > '(1x,e16.8))'
823  write(lunlist, fmt)
824  + context, rname, chset, chdet, nh, (chnam(k), k=1,nh),
825  + (nbits(k), k=1,nh), (orig(k), k=1,nh), (fact(k), k=1,nh)
826  endif
827  if (luncode.ne.0) then
828  write(luncode,'(''{'')')
829  call astocp('chnam',chnam,nh)
830  call aitocp('nbits',nbits,nh)
831  call artocp('orig',orig,nh)
832  call artocp('fact',fact,nh)
833  write(luncode,1000) chset,chdet,nh
834  1000 format('G4gsdeth(chset="',a,'",chdet="',a,'",nh=',i4,
835  + ',chnam,nbits,orig,fact);')
836  write(luncode,'(''}'')')
837  endif
838 *
839  end
840 *
841 #ifndef CALL_GEANT
842  subroutine gsdetd(chset, chdet, nd, chnam, nbits)
843 #else
844  subroutine ksdetd(chset, chdet, nd, chnam, nbits)
845 #endif
846 ************************************************************************
847 ************************************************************************
848  implicit none
849  integer nd, nbits(nd), k
850  character rname*6, chset*4, chdet*4, chnam(nd)*4, fmt*150
851 #include "G3toG4.inc"
852  data rname /'GSDETD'/
853 *
854  call check_lines
855 #ifdef CALL_GEANT
856  if (dogeom) call gsdetd(chset, chdet, nd, chnam, nbits)
857 #endif
858  if (lunlist.ne.0) then
859 * write(lunlist,
860 * + '(a4,1x,a6,1x,a4,1x,a4,i5,<nd>(1x,a4),<nd>i5)')
861 * + context, rname, chset, chdet, nd, (chnam(k), k=1,nd),
862 * + (nbits(k), k=1,nd)
863  write(fmt,'(A,I3,A,I3,A)')
864  + '(a4,1x,a6,1x,a4,1x,a4,i5,',max(nd,1),'(1x,a4),',
865  > max(nd,1),'i5)'
866  write(lunlist,fmt)
867  + context, rname, chset, chdet, nd, (chnam(k), k=1,nd),
868  + (nbits(k), k=1,nd)
869  endif
870  if (luncode.ne.0) then
871  write(luncode,'(''{'')')
872  call astocp('chnam',chnam,nd)
873  call aitocp('nbits',nbits,nd)
874  write(luncode,1000) chset, chdet, nd
875  1000 format('G4gsdetd(chset="',a,'",chdet="',a,'",nd=',i4,
876  + ',chnam,nbits);')
877  write(luncode,'(''}'')')
878  endif
879 *
880  end
881 *
882 #ifndef CALL_GEANT
883  subroutine gsdetu(chset, chdet, nupar, upar)
884 #else
885  subroutine ksdetu(chset, chdet, nupar, upar)
886 #endif
887 ************************************************************************
888 ************************************************************************
889  implicit none
890  integer nupar, k
891  real upar(nupar)
892  character rname*6, chset*4, chdet*4, fmt*150
893 #include "G3toG4.inc"
894  data rname /'GSDETU'/
895 *
896  call check_lines
897 #ifdef CALL_GEANT
898  if (dogeom) call gsdetu(chset, chdet, nupar, upar)
899 #endif
900  if (lunlist.ne.0) then
901 * write(lunlist,
902 * + '(a4,1x,a6,1x,a4,1x,a4,i5,<nupar>e15.8)')
903 * + context, rname, chset, chdet, nupar, (upar(k), k=1,nupar)
904  write(fmt,'(A,I3,A)')
905  + '(a4,1x,a6,1x,a4,1x,a4,i5,',max(nupar,1),'(1x,e16.8))'
906  write(lunlist,fmt)
907  + context, rname, chset, chdet, nupar, (upar(k), k=1,nupar)
908  endif
909  if (luncode.ne.0) then
910  write(luncode,'(''{'')')
911  call g3ldpar(upar,nupar)
912  write(luncode,1000) chset, chdet, nupar
913  1000 format('G4gsdetu(chset="',a,'",chdet="',a,'",npar=',
914  + i4,',par);')
915  write(luncode,'(''}'')')
916  endif
917 *
918  end
919 *
920 #ifndef CALL_GEANT
921  subroutine ggclos
922 #else
923  subroutine kgclos
924 #endif
925 ************************************************************************
926 ************************************************************************
927  implicit none
928  character rname*6
929 #include "G3toG4.inc"
930  data rname /'GGCLOS'/
931 *
932  call check_lines
933 #ifdef CALL_GEANT
934  if (dogeom) call ggclos
935 #endif
936  if (lunlist.ne.0) then
937  write(lunlist,'(a4,1x,a6)') context, rname
938  close(lunlist)
939  endif
940  if (luncode.ne.0) then
941  write(luncode,'(''//GeoMgr->CloseGeometry();'')')
942  write(luncode,'(''}'')')
943  call g3main
944  close(luncode)
945  endif
946 *
947  end
948 
949  subroutine checkshape(name, shape, par, npar)
950  implicit none
951 ************************************************************************
952 * convert TRAP, PARA and GTRA to external form
953 ************************************************************************
954  character name*4, shape*4
955  real ph, par(*), tt, raddeg
956  integer npar
957 
958  raddeg = 180./3.1415926
959 
960  if (shape(1:3).eq.'BOX'.and.npar.ne.3) then
961  print *,'!! error, BOX with ',npar,' parameters, vol ',name
962  endif
963  if (shape.eq.'TRD1'.and.npar.ne.4) then
964  print *,'!! error, TRD1 with ',npar,' parameters, vol ',name
965  endif
966  if (shape.eq.'TRD2'.and.npar.ne.5) then
967  print *,'!! error, TRD2 with ',npar,' parameters, vol ',name
968  endif
969  if (shape.eq.'TRAP'.and.npar.ne.35.and.npar.ne.11) then
970 *** G3 sets 11 to 35. Why?
971  print *,'!! error, TRAP with ',npar,' parameters, vol ',name
972  endif
973  if (shape.eq.'TUBE'.and.npar.ne.3) then
974  print *,'!! error, TUBE with ',npar,' parameters, vol ',name
975  endif
976  if (shape.eq.'TUBS'.and.npar.ne.5) then
977  print *,'!! error, TUBS with ',npar,' parameters, vol ',name
978  endif
979  if (shape.eq.'CONE'.and.npar.ne.5) then
980  print *,'!! error, CONE with ',npar,' parameters, vol ',name
981  endif
982  if (shape.eq.'CONS'.and.npar.ne.7) then
983  print *,'!! error, CONS with ',npar,' parameters, vol ',name
984  endif
985  if (shape.eq.'SPHE'.and.npar.ne.6) then
986  print *,'!! error, SPHE with ',npar,' parameters, vol ',name
987  endif
988  if (shape.eq.'PARA'.and.npar.ne.6) then
989  print *,'!! error, PARA with ',npar,' parameters, vol ',name
990  endif
991  if (shape.eq.'PARA') then
992 *
993 * ** PARA
994 *
995  ph = 0.
996  if (par(5).ne.0.) ph = atan2(par(6),par(5))*raddeg
997  tt = sqrt(par(5)**2+par(6)**2)
998  par(4) = atan(par(4))*raddeg
999  if (par(4).gt.90.0) par(4) = par(4)-180.0
1000  par(5) = atan(tt)*raddeg
1001  if (ph.lt.0.0) ph = ph + 360.0
1002  par(6) = ph
1003  end if
1004  if (shape.eq.'TRAP') then
1005 *
1006 * ** TRAP
1007 *
1008  npar=11
1009  ph = 0.
1010  if (par(2).ne.0.) ph = atan2(par(3),par(2))*raddeg
1011  tt = sqrt(par(2)**2+par(3)**2)
1012  par(2) = atan(tt)*raddeg
1013  if (ph.lt.0.0) ph = ph+360.0
1014  par(3) = ph
1015  par(7) = atan(par(7))*raddeg
1016  if (par(7).gt.90.0) par(7) = par(7)-180.0
1017  par(11)= atan(par(11))*raddeg
1018  if (par(11).gt.90.0) par(11) = par(11)-180.0
1019 
1020  end if
1021  end