14 c call hbook1(nh0+1,' ngam$',ngmax,0.,fngmax,0)
15 c call hbook1(nh0+2,' egam$',100,0.,egmax,0)
16 c call hbook1(nh0+3,' etot$',150,0.,etmax,0)
17 c call hbook1(nh0+4,' pt $',100,0.,ptmax,0)
18 c call hbook2(nh0+5,' x*y $',50,-smax,smax,50,-smax,smax,0)
19 c call hbook1(nh0+6,' x $',500,-smax,smax,0)
20 c call hbook1(nh0+7,' y $',500,-smax,smax,0)
21 c call hbook1(nh0+8,' ptalice $',100,0.,ptmax,0)
25 c call hf1(nh0+1,float(kg),1.)
29 pt=sqrt(p4(1,k)**2+p4(2,k)**2)
30 ptalice=sqrt(p4(1,k)**2+p4(3,k)**2)
31 c call hf1(nh0+2,e(k),1.)
32 c call hf2(nh0+5,xw(k),yw(k),1.)
33 c call hf1(nh0+6,xw(k),1.)
34 c call hf1(nh0+7,yw(k),1.)
35 c call hf1(nh0+4,pt,1.)
36 c call hf1(nh0+8,ptalice,1.)
38 c call hf1(nh0+3,etot,1.)
44 common /comgeom/igeomflag
47 include 'comwlgen.for'
48 external xgamvert,ygamvert,zgamvert
56 if(igeomflag.eq.2) then
57 rg=zgdev(k)+xyzwall(3,igd)
66 c write (*,*) ' k zgdev xyzwall zv',k,zgdev(k),xyzwall(3,igd),zv
67 zg=zgdev(k)+xyzwall(3,igd)-zv
69 rr=sqrt(xg**2+yg**2+zg**2)
79 real function xgamvert(ivert)
83 elseif(ivert.gt.ntotvertx) then
85 write (*,*) ' wrong vertex number '
87 xgamvert=gamvertex(1,ivert)
92 real function ygamvert(ivert)
96 elseif(ivert.gt.ntotvertx) then
98 write (*,*) ' wrong vertex number '
100 ygamvert=gamvertex(2,ivert)
105 real function zgamvert(ivert)
109 elseif(ivert.gt.ntotvertx) then
111 write (*,*) ' wrong vertex number '
113 zgamvert=gamvertex(3,ivert)
118 subroutine clrcompart
119 include 'compart.for'
124 subroutine addp4tocompart(ncrad)
126 common /comgeom/igeomflag
127 include 'event_format.inc'
129 include 'compart.for'
138 nmmmax=crystals_matrix_amount_PHOS
141 fir(i)=fi0*float(nmmmax-1)-2.*fi0*float(i-1)
147 if(nbpart+1.gt.nbpartmax) then
148 write (*,*) ' incr nbpartmax ',nbpartmax
153 if(igeomflag.eq.1) then
154 call rotrefr(p4(1,k),ppart(1,nbpart),fir(ncrad))
156 call rotrefc(p4(1,k),ppart(1,nbpart))
163 subroutine extrp4wrf(ncrad)
165 common /comgeom/igeomflag
166 include 'event_format.inc'
168 include 'compart.for'
169 include 'comggen.for'
181 nmmmax=crystals_matrix_amount_PHOS
183 fir(i)=fi0*float(nmmmax-1)-2.*fi0*float(i-1)
189 c if(igeomflag.eq.1) then
190 c call rotrefr_a(p4temp,ppart(1,n),fir(ncrad))
192 c call rotrefc_a(p4temp,ppart(1,n))
196 c call xydetwlnew(xx,yy,p4temp)
197 c call mdetwlnew(mmww,xx,yy,ijdd,1)
202 if(igeomflag.eq.1) then
203 call rotrefr_a(p4temp,g(1,n),fir(ncrad))
205 call rotrefc_a(p4temp,g(1,n))
209 call xydetwlnew(xx,yy,p4temp)
210 write (*,*) ' xx,yy ',xx,yy,xgamgen(n),ygamgen(n)
211 c call mdetwlnew(mmww,xx,yy,ijdd,1)
218 subroutine calcgamgennew(ncrad)
220 common /comgeom/igeomflag
221 include 'event_format.inc'
223 include 'compart.for'
224 include 'comggen.for'
236 nmmmax=crystals_matrix_amount_PHOS
238 fir(i)=fi0*float(nmmmax-1)-2.*fi0*float(i-1)
244 call ucopy(g(1,n),p4ggen(1,n),3)
245 p4ggen(4,n)=sqrt(g(1,n)**2+g(2,n)**2+g(3,n)**2)
246 egamgen(n)=p4ggen(4,n)
247 if(igeomflag.eq.1) then
248 call rotrefr_a(p4temp,p4ggen(1,n),fir(ncrad))
250 call rotrefc_a(p4temp,p4ggen(1,n))
254 call xydetwlnew(xx,yy,p4temp)
267 rx=xgamgen(n1)-xgamgen(n2)
268 ry=ygamgen(n1)-ygamgen(n2)
270 if(rr.lt.rmingen(n1)) then
274 if(rr.lt.rmingen(n2)) then
285 subroutine xydetwlnew(x,y,p4)
286 common /comgeom/igeomflag
291 if(igeomflag.eq.2) then
294 fi1=atan(p4(1)/p4(3))
295 g12=sqrt(p4(3)**2+p4(1)**2)
298 elseif(igeomflag.eq.1) then
303 write (*,*) ' check geometry flag '
311 subroutine rotrefr(p1,p2,fi)
318 p2(1)=xp*cosfi+yp*sinfi
319 p2(2)=yp*cosfi-xp*sinfi
325 subroutine rotrefc(p1,p2)
334 subroutine rotrefr_a(p1,p2,fi)
340 p1(3)=p2(1)*cosfi-p2(2)*sinfi
341 p1(1)=p2(2)*cosfi+p2(1)*sinfi
348 subroutine rotrefc_a(p1,p2)