8d0a87de7c860f50cd3790569b84d17ee4331c6c
[u/mrichter/AliRoot.git] / EPOS / epos167 / epos-tim-155.f
1 c---------------------------------------------------------------------
2       subroutine timann   
3 c---------------------------------------------------------------------
4 c      electron-positron
5 c---------------------------------------------------------------------
6       include 'epos.inc'
7       include 'epos.incsem'
8       common/nfla/nfla 
9       parameter (ntim=1000)
10       common/cprt/nprtj,pprt(5,ntim),idprt(ntim),iorprt(ntim)
11      &,idaprt(2,ntim)
12       integer jcp(6,2),jcm(6,2)
13       dimension jorprt(ntim),p1(5)
14
15       call utpri('timann',ish,ishini,5)
16       egyevt=engy
17       qsqevt=engy**2
18
19  123  nprtj=1
20       en=engy
21       q20=0.
22       nptl=0
23       nptl=nptl+1               !electron
24       pptl(1,nptl)=0.
25       pptl(2,nptl)=0.
26       pptl(3,nptl)=sqrt(en**2/4-2.61121e-7)
27       pptl(4,nptl)=en/2.
28       pptl(5,nptl)=0.000511
29       idptl(nptl)=12
30       istptl(nptl)=1
31
32       nptl=nptl+1               !positron
33       pptl(1,nptl)=0.
34       pptl(2,nptl)=0.
35       pptl(3,nptl)=-sqrt(en**2/4-2.61121e-7)
36       pptl(4,nptl)=en/2.
37       pptl(5,nptl)=0.000511
38       idptl(nptl)=-12
39       istptl(nptl)=1
40
41       nptl=nptl+1               !virtual gamma
42       pptl(1,nptl)=0.
43       pptl(2,nptl)=0.
44       pptl(3,nptl)=0.
45       pptl(4,nptl)=en
46       pptl(5,nptl)=en
47       idptl(nptl)=10
48       iorptl(nptl)=nptl-2
49       jorptl(nptl)=nptl-1
50       istptl(nptl)=1
51
52       pprt(1,nprtj)=0.
53       pprt(2,nprtj)=0.
54       pprt(3,nprtj)=0.
55       pprt(4,nprtj)=en    
56       pprt(5,nprtj)=en    
57       idaprt(1,nprtj)=2
58       idaprt(2,nprtj)=3
59       if(q20.gt.0.)then
60         pprt(5,nprtj)=sqrt(q20)
61       else
62         pprt(5,nprtj)=en
63       endif
64       nfla=0
65       do i=1,naflav
66         call idmass(i,am)
67         if (2.*am.lt.en) nfla=i
68       enddo
69       if(itflav.eq.0)then
70         s=engy**2 
71         dlz=2.4
72         amz=91.1885
73         al=1./real(137.035989d0)
74         gf=1.16639e-5
75         ak=sqrt(2.)*gf*amz**2/(16*pi*al)
76         chi1=ak*s*(s-amz**2)/((s-amz**2)**2+dlz**2*amz**2)
77         chi2=ak**2*s**2/((s-amz**2)**2+dlz**2*amz**2)
78         qe=-1.
79         ve=0.25-2.*qe*0.232
80         ae=-0.5
81         qf=2./3.
82         vf=sign(.5,qf)-2.*qf*0.232 
83         af=sign(.5,qf)
84         dsmax1=
85      $       2.*(qf**2-2.*qf*ve*vf*chi1 
86      $       +(ae**2+ve**2)*(af**2+vf**2)*chi2)
87      $       + abs(-4*qf*ae*af*chi1+8*ae*ve*af*vf*chi2)
88     
89         qf=-1./3.
90         vf=sign(.5,qf)-2.*qf*0.232 
91         af=sign(.5,qf)
92         dsmax2=
93      $       2.*(qf**2-2.*qf*ve*vf*chi1 
94      $       +(ae**2+ve**2)*(af**2+vf**2)*chi2)
95      $       + abs(-4*qf*ae*af*chi1+8*ae*ve*af*vf*chi2)
96
97  100    iq1=1+INT(nfla*rangen()) 
98         call idchrg(iq1,qf)
99         ct=-1.+2.*rangen()
100         vf=sign(.5,qf)-2.*qf*0.232 
101         af=sign(.5,qf)
102
103         dsigma=
104      $       (1.+ct**2)*(qf**2-2.*qf*ve*vf*chi1 
105      $       +(ae**2+ve**2)*(af**2+vf**2)*chi2)
106      $       + ct*(-4*qf*ae*af*chi1+8*ae*ve*af*vf*chi2)
107         if(rangen().gt.dsigma/max(dsmax1,dsmax2)) goto 100
108       else
109         iq1=itflav
110       endif
111       if(rangen().lt.0.5)iq1=-iq1
112       
113       nprtj=nprtj+1
114       idprt(nprtj)=iq1
115       pprt(4,nprtj)=en/2.
116       pprt(5,nprtj)=en
117       iorprt(nprtj)=1
118       
119       nprtj=nprtj+1
120       idprt(nprtj)=-iq1
121       pprt(4,nprtj)=en/2.
122       pprt(5,nprtj)=en
123       iorprt(nprtj)=1
124
125       call timsho(2,3)
126
127       jorprt(1)=0               !!color-connection, no origin!!
128       jt=1
129       if(idprt(idaprt(1,1)).lt.0)jt=2
130       do i=1,nprtj 
131         if(idaprt(1,i).ne.0) then
132           js=jt
133           if(idprt(i).lt.0.and.
134      &      ((idprt(idaprt(2,i)).eq.9.and.jt.eq.1).or.
135      &      (idprt(idaprt(1,i)).eq.9.and.jt.eq.2)))then 
136             js=3-jt
137           elseif(idprt(i).gt.0.and.idprt(i).ne.9.and.
138      &        ((idprt(idaprt(2,i)).eq.9.and.jt.eq.2).or.
139      &        (idprt(idaprt(1,i)).eq.9.and.jt.eq.1)))then 
140             js=3-jt
141           elseif(idprt(i).eq.9.and.idprt(idaprt(1,i)).ne.9.and.
142      &        ((idprt(idaprt(1,i)).lt.0.and.jt.eq.2).or.
143      &        (idprt(idaprt(1,i)).gt.0.and.jt.eq.1)))then
144             js=3-jt
145           endif
146           jorprt(idaprt(3-js,i))=jorprt(i)
147           jorprt(i)=idaprt(js,i)
148           jorprt(idaprt(js,i))=idaprt(3-js,i)
149         endif
150       enddo
151       
152       if(ish.ge.5)then
153         i=1    
154         do while(i.ne.0)
155           if(idaprt(1,i) .eq. 0 ) then 
156             write(ifch,*) idprt(i)
157             write(ifch,*) '|'
158           endif
159           i=jorprt(i)
160         enddo
161       endif
162
163       iptl=nptl
164       i=1    
165       do while(i.gt.0)
166         if(idaprt(1,i) .eq. 0) then 
167           nptl=nptl+1
168           do j=1,5
169             pptl(j,nptl)=pprt(j,i)
170           enddo
171           idptl(nptl)=idprt(i)
172           istptl(nptl)=20
173           iorptl(nptl)=iptl
174           qsqptl(nptl)=en**2
175         endif
176         i=jorprt(i)
177       enddo
178       ifrptl(1,iptl)=iptl+1
179       ifrptl(2,iptl)=nptl
180
181       nk1=iptl+1
182  441  nk2=nk1+1
183       do while (idptl(nk2).eq.9)
184         nk2=nk2+1
185       enddo
186       do i=1,4
187         p1(i)=0.
188         do j=nk1,nk2
189           p1(i)=p1(i)+pptl(i,j)
190         enddo
191       enddo
192       p1(5)=sqrt(max(0.,p1(4)**2-p1(3)**2-p1(2)**2-p1(1)**2))
193       do i=1,2
194         do j=1,6
195           jcm(j,i)=0
196           jcp(j,i)=0
197         enddo
198       enddo
199       ii=1
200       if(idptl(nk1).lt.0)ii=2
201       jcp(abs(idptl(nk1)),ii)=1
202       jcm(abs(idptl(nk2)),3-ii)=1
203       amm= utamnx(jcp,jcm)
204       if(amm.gt.p1(5))goto 123
205       nk1=nk2+1
206       if(nk1.lt.nptl)goto 441
207
208
209       if(ish.ge.5)then
210         write(ifch,*)
211         do i=1,nprtj
212           write(ifch,98) i,(pprt(j,i),j=1,5),idprt(i)
213      &    ,iorprt(i),idaprt(1,i),idaprt(2,i),jorprt(i)
214         enddo
215         write(ifch,*)
216         do i=1,nprtj
217           if(pprt(5,i).eq.0.) 
218      &    write(ifch,99) i,(pprt(j,i),j=1,5),idprt(i)
219         enddo
220         write(ifch,*)
221         write(ifch,*)
222       endif
223       
224  99   format(i4,5g10.3,1i4)
225  98   format(i4,5g10.3,5i4)
226             
227       call utprix('timann',ish,ishini,5)
228       return
229       end
230       
231 c---------------------------------------------------------------------
232       subroutine timsh1(q20,en,idfla)
233 c---------------------------------------------------------------------
234       include 'epos.inc'
235       include 'epos.incsem'
236       common/nfla/nfla
237       parameter (ntim=1000)
238       common/cprt/nprtj,pprt(5,ntim),idprt(ntim),iorprt(ntim)
239      &,idaprt(2,ntim)
240       common/cprtx/nprtjx,pprtx(5,2)
241       nfla=naflav
242       nprtj=1
243       pprt(1,nprtj)=0.
244       pprt(2,nprtj)=0.
245       pprt(3,nprtj)=0.
246       pprt(4,nprtj)=en    
247       pprt(5,nprtj)=sqrt(q20)
248       idprt(nprtj)=idfla
249       call timsho(1,0)
250       nprtjx=0
251       return
252       end
253
254
255 c---------------------------------------------------------------------
256       subroutine timsh2(q20,q21,en,idfla1,idfla2) 
257 c---------------------------------------------------------------------
258       include 'epos.inc'
259       include 'epos.incsem'
260       common/nfla/nfla
261       parameter (ntim=1000)
262       common/cprt/nprtj,pprt(5,ntim),idprt(ntim),iorprt(ntim)
263      &,idaprt(2,ntim)
264       common/cprtx/nprtjx,pprtx(5,2)
265       nfla=naflav
266
267       nprtj=1
268       pprt(1,nprtj)=0.
269       pprt(2,nprtj)=0.
270       pprt(3,nprtj)=0.
271       pprt(4,nprtj)=en/2. 
272       pprt(5,nprtj)=sqrt(q20)
273       idprt(nprtj)=idfla1
274
275       nprtj=2
276       pprt(1,nprtj)=0.
277       pprt(2,nprtj)=0.
278       pprt(3,nprtj)=0.
279       pprt(4,nprtj)=en/2.
280       pprt(5,nprtj)=sqrt(q21)
281       idprt(nprtj)=idfla2
282
283       call timsho(1,2)
284
285       nprtjx=1
286       pprtx(1,nprtjx)=0.
287       pprtx(2,nprtjx)=0.
288       pprtx(3,nprtjx)=en/2. 
289       pprtx(4,nprtjx)=en/2. 
290       pprtx(5,nprtjx)=0
291       nprtjx=2
292       pprtx(1,nprtjx)=0.
293       pprtx(2,nprtjx)=0.
294       pprtx(3,nprtjx)=-en/2. 
295       pprtx(4,nprtjx)=en/2.
296       pprtx(5,nprtjx)=0
297
298       return
299       end
300
301 c---------------------------------------------------------------------
302       subroutine timsho(j1,j2)
303 c---------------------------------------------------------------------
304
305
306       include 'epos.inc'
307       include 'epos.incsem'
308       parameter (ntim=1000)
309       common/cprt/nprtj,pprt(5,ntim),idprt(ntim),iorprt(ntim)
310      &,idaprt(2,ntim)
311       dimension pz(ntim),id(2,ntim)
312       dimension ij(2),ee(2),amm2(-6:10)
313       logical first
314       common/nfla/nfla
315       call utpri('timsho',ish,ishini,5)
316       if(iappl.eq.6)then
317         do i=1,6
318           call idmass(i,am)
319           amm2(i)=am**2
320           amm2(-i)=am**2
321         enddo
322       else
323         do i=1,6
324           amm2(i)=0.
325           amm2(-i)=0.
326         enddo
327       endif
328       amm2(9)=0.
329       amm2(0)=0.
330       amm2(10)=0.
331       nn=nprtj
332       ij(1)=j1
333       ij(2)=j2
334       ii2=2
335       if(j2.eq.0)ii2=1
336
337       ii=1
338       first=.true.
339       n10=0
340  10   n10=n10+1
341       if(float(n10).gt.1e7)then
342         goto9999
343       endif
344       io=iorprt(ij(ii))
345       idfl=idprt(ij(ii))
346       q2start=pprt(5,ij(ii))**2
347       E=pprt(4,ij(ii))
348       if(ij(2).eq.j2.and.ii2.eq.2)E=pprt(4,ij(1))+pprt(4,ij(2))
349       zetamx=0.
350       if(ij(1).ne.j1)then
351         zetamx=pprt(5,io)/pprt(4,io)/sqrt(pz(io)*(1.-pz(io))) 
352       endif
353
354 c      call timdev(idfl,q2start,E,zetamx,idfla,idflb,qa2,z)
355 c......................................................................
356  380  q2=q2start
357       z=0.5
358       PT2MIN=max(qcdlam*1.1,q2fin)
359       ALFM=LOG(PT2MIN/qcdlam) 
360       if(ish.ge.9)then
361         write (ifch,*) '---------------------',ii
362      $       ,pprt(5,ij(1))     !,pprt(5,ij(2))
363         write(ifch,*) ' idfl,q2start,zetamx:',idfl,q2start,zetamx
364       endif
365       if (q2.lt.4.*q2fin+2.*amm2(idfl) )then
366         if(ish.ge.9)then
367           write(ifch,'(a,i4,i4,2f15.8)') 'null:',0.
368         endif
369         q2=amm2(idfl)
370         idfla=0
371         idflb=0
372         goto 999
373       endif
374       
375  390  zc=.5*(1.-sqrt(max(0.000001,1.-4.*q2fin/q2)))
376       if(ish.ge.9)then
377         write(ifch,*) 'zc=',zc
378       endif
379       IF(idfl.EQ.9) THEN 
380         FBR=6.*LOG((1.-ZC)/ZC)+nfla*(0.5-ZC) 
381       ELSE
382         FBR=(8./3.)*LOG((1.-ZC)/ZC) 
383       endif
384       
385       B0=(33.-2.*nfla)/6.
386       
387       
388 c.....select new q2
389       r=rangen()
390       q2=q2*exp(log(r)*B0*ALFM/FBR)
391       if(ish.ge.9)then
392         write(ifch,*) 'q^2=',q2
393       endif
394       if (q2.lt.4.*q2fin+2.*amm2(idfl))then
395         q2=amm2(idfl)
396         if(ish.ge.9)then
397           write(ifch,'(a,i4,i4,2f15.8)') 'null:',0.
398         endif
399         idfla=0
400         idflb=0
401         goto 999
402       endif
403       
404 c.....select flavor and z-value .....................................
405       IF(idfl.EQ.9) THEN 
406         if(rangen()*FBR.lt.nfla*(0.5-ZC))then
407                                 ! .................g -> qqbar
408           Z=ZC+(1.-2.*ZC)*rangen()
409           IF(Z**2+(1.-Z)**2.LT.rangen()) GOTO 390 
410           idfla=int(1.+rangen()*real(nfla))
411           idflb=-idfla
412         else                    !..................g -> gg
413           Z=(1.-ZC)*(ZC/(1.-ZC))**rangen()
414           IF(rangen().GT.0.5) Z=1.-Z 
415           IF((1.-Z*(1.-Z))**2.lt.rangen()) GOTO 390 
416           idfla=9
417           idflb=9
418         endif
419       ELSE
420         Z=1.-(1.-ZC)*(ZC/(1.-ZC))**rangen() !!........q -> qg
421         IF(1.+Z**2.LT.2.*rangen()) GOTO 390 
422           idfla=idfl
423           idflb=9
424       endif
425
426       if(alfm/log(q2*z*(1.-z)/qcdlam).lt.rangen()) goto 390 
427       
428       if(ij(1).ne.j1.or.ij(2).eq.0)then
429         if(E.lt.sqrt(q2))goto 390
430         pzz=sqrt((E-sqrt(q2))*(E+sqrt(q2)))
431         pt2=(E**2*(z*(1.-z)*q2-z*amm2(idflb)-(1.-z)*amm2(idfla))
432      $       -.25*(amm2(idflb)-amm2(idfla)-q2)**2+q2*amm2(idfla))/pzz**2
433         if (pt2.le.0.) then
434           if(ish.ge.9)then
435             write(ifch,*) 'z not good for pt2:',z,pt2
436           endif
437           goto 390
438         endif
439       endif
440
441       zeta = sqrt(q2)/E/sqrt(z*(1.-z))
442       if(zetamx.gt.0)then
443         if (zeta.gt.zetamx)then
444           if(ish.ge.9)then
445             write(ifch,*) zeta,' > ',zetamx,'zeta-Ablehnung'
446           endif
447           goto 390
448         endif
449       endif
450  999  continue
451 c......................................................................
452     
453       pprt(5,ij(ii))=sqrt(q2)
454       id(1,ii)=idfla
455       id(2,ii)=idflb
456       pz(ij(ii))=z
457       if(first)then
458         ii=2
459         first=.false.
460         if(ii2.eq.2)goto 10
461       endif
462
463       if(ij(2).eq.0)then
464         z1=1.
465         pt=0.
466         alpha=0.
467         pprt(3,ij(1))=sqrt(E**2-pprt(5,ij(1))**2)
468       elseif(ij(1).eq.j1.and.ij(2).eq.j2)then
469         E=pprt(4,ij(1))+pprt(4,ij(2))
470         ee(1)=E*.5+(pprt(5,ij(1))**2-pprt(5,ij(2))**2)/2./E
471         ee(2)=E-ee(1)
472         ii=1
473         do while (ii.le.2)
474           if(ee(ii)-pprt(5,ij(ii)).lt.0.) then
475             if(ish.ge.5)then
476               write(ifch,*) 'goto 11'
477             endif
478             ii=1
479             if ( pprt(5,ij(1))**2-amm2(idprt(ij(1))).lt.
480      $           pprt(5,ij(2))**2-amm2(idprt(ij(2))) ) ii=2
481             goto 10
482           endif
483 c          zc=.5*(1.-sqrt(1.-pprt(5,ij(ii))**2/ee(ii)**2))
484 c          if(pz(ij(ii)).lt.zc.or.pz(ij(ii)).gt.1.-zc)then
485 c            if(ish.ge.7)then
486 c              write(ifch,*) 'first branching rejected'
487 c            endif
488 c            goto 10
489 c          endif
490           z=pz(ij(ii))
491           q2=pprt(5,ij(ii))**2
492           pzz=sqrt((ee(ii)-pprt(5,ij(ii)))*(ee(ii)+pprt(5,ij(ii))))
493           if(pzz.gt.0.)then
494            pt2=(ee(ii)**2*(z*(1.-z)*q2-z*amm2(id(2,ii))
495      $          -(1.-z)*amm2(id(1,ii)))
496      $          -.25*(amm2(id(2,ii))-amm2(id(1,ii))-q2)**2
497      $          +q2*amm2(id(1,ii)))/pzz**2
498           else
499            pt2=0.
500           endif
501           if(id(1,ii).ne.0.and.pt2.le.0.)then
502             if(ish.ge.7)then
503               write(ifch,*) 'first branching rejected for pt2',ii
504      $             ,z1,q2,ee(ii),pprt(5,ij(ii)),id(1,ii),id(2,ii)
505             endif
506             goto 10
507           endif
508           ii=ii+1
509         enddo
510         z1=ee(1)/E
511         z2=ee(2)/E
512         if(ish.ge.7)then
513           write(ifch,*) 'z of first branching',z1
514         endif
515         pprt(3,ij(1))= sqrt(max(0.,ee(1)**2-pprt(5,ij(1))**2))
516         pprt(3,ij(2))=-sqrt(max(0.,ee(2)**2-pprt(5,ij(2))**2))
517         pt=0.
518         alpha=0.
519       else
520         E=pprt(4,io)
521         z1=pz(io)
522         z2=1.-z1
523         am0=pprt(5,io)
524         am1=pprt(5,ij(1))
525         am2=pprt(5,ij(2))
526         aM=am2**2-am0**2-am1**2
527         pzz=sqrt((E-am0)*(E+am0))
528         pprt(3,ij(1))=.5*(aM+2.*z1*E**2)/pzz
529         pprt(3,ij(2))=pzz-pprt(3,ij(1))
530         pt2=(E**2*(z1*z2*am0**2-z1*am2**2-z2*am1**2)
531      $       -.25*aM**2+am0**2*am1**2)/pzz**2
532         if(ish.ge.8)then
533           write(ifch,*) 'pt2,pzz=',pt2,pzz,z1**2*E*pprt(5,io),z1,E
534      $         ,pprt(5,io)
535         endif
536  111    if(pt2.lt.0.) then 
537           ii=1
538           if ( pprt(5,ij(1))**2-amm2(idprt(ij(1))).lt.
539      $         pprt(5,ij(2))**2-amm2(idprt(ij(2))) ) ii=2
540           goto 10 
541         endif
542         pt=sqrt(pt2)
543         alpha=2.*pi*rangen()
544       endif
545       pprt(4,ij(1))=z1*E
546       pprt(1,ij(1))=cos(alpha)*pt
547       pprt(2,ij(1))=sin(alpha)*pt
548       if(ii2.eq.2)then
549         pprt(4,ij(2))=z2*E
550         pprt(1,ij(2))=-cos(alpha)*pt
551         pprt(2,ij(2))=-sin(alpha)*pt
552       endif
553       if(ij(1).ne.j1.and.(pprt(1,io).ne.0..or.pprt(2,io).ne.0.))then
554         do ii=1,2
555           call utrota(-1,pprt(1,io),pprt(2,io),pprt(3,io)
556      &    ,pprt(1,ij(ii)),pprt(2,ij(ii)),pprt(3,ij(ii)))
557         enddo
558       endif
559       if(ij(1).ne.j1.and.pprt(3,io).lt.0.)then
560         do k=1,3
561           do ii=1,2
562             pprt(k,ij(ii)) = -pprt(k,ij(ii))
563           enddo
564         enddo
565       endif
566       do ii=1,ii2
567         if(id(1,ii).ne.0)then
568           idprt(nprtj+1)=id(1,ii)
569           idprt(nprtj+2)=id(2,ii)
570           pprt(4,nprtj+1)=pz(ij(ii))*pprt(4,ij(ii))
571           pprt(5,nprtj+1)=pz(ij(ii))*pprt(5,ij(ii))
572           pprt(4,nprtj+2)=(1.-pz(ij(ii)))*pprt(4,ij(ii))
573           pprt(5,nprtj+2)=(1.-pz(ij(ii)))*pprt(5,ij(ii))
574           iorprt(nprtj+1)=ij(ii)
575           iorprt(nprtj+2)=ij(ii)
576           idaprt(1,ij(ii))=nprtj+1
577           idaprt(2,ij(ii))=nprtj+2
578           idaprt(1,nprtj+1)=0
579           idaprt(2,nprtj+1)=0
580           idaprt(1,nprtj+2)=0
581           idaprt(2,nprtj+2)=0
582           nprtj=nprtj+2
583         else
584           idaprt(1,ij(ii))=0
585           idaprt(2,ij(ii))=0
586         endif
587       enddo
588       if(ish.ge.5)then
589         if(ij(1).ne.j1)then
590           write(ifch,98) io,(pprt(j,io),j=1,5),pz(io)
591      &    ,idprt(io)
592           write(ifch,*) '->'
593         endif
594         write(ifch,99) ij(1),(pprt(j,ij(1)),j=1,5),pz(ij(1))
595      &  ,idprt(ij(1)),'->',id(1,1),id(2,1)
596         if(ij(2).ne.0)write(ifch,99) ij(2),
597      &  (pprt(j,ij(2)),j=1,5),pz(ij(2))
598      &  ,idprt(ij(2)),'->',id(1,2),id(2,2)
599         write(ifch,*)
600       endif
601
602  98   format(i4,6g10.3,1i4)
603  99   format(i4,6g10.3,1i4,a,2i4)
604
605
606       if(ij(1).le.nn)ij(1)=nn-1
607       ij(1)=ij(1)+2
608       ij(2)=ij(1)+1
609       ii=1
610       ii2=2
611       first=.true.
612       if(ij(1).le.nprtj)goto10
613
614       if(ish.ge.5)then
615         write(ifch,*)
616         do i=1,nprtj
617           write(ifch,'(i4,1x,5g10.4,2i4,a,2i4)') 
618      &    i,(pprt(j,i),j=1,5),idprt(i)
619      &    ,iorprt(i),'  -->',idaprt(1,i),idaprt(2,i)
620         enddo
621         write(ifch,*)
622       endif
623  9999 call utprix('timsho',ish,ishini,5)
624       return
625       end
626