program NOTEND9 ! ! Removes good tracks from ICOOL file FOR009.DAT ! Writes file containing all events that don't make it to the end ! parameter ( mxtrk = 100000 ) real*8 tp,xp(3),pp(3),bfld(3),evtwt,efld(3),sarc,pol(3) integer igd(mxtrk),ibd(mxtrk) character*80 line,fmt data igd/mxtrk*0/,ibd/mxtrk*0/,idec/0/,iexp/0/ ! version = 1.06 write(*,'(a,f6.2)') 'Welcome to program NOTEND9, version ',version ! open(unit=9,file='for009.dat',status='old',iostat=ioc) if(ioc .ne. 0) then write(*,*)'Couldnt open FOR009.DAT' go to 900 end if open(unit=10,file='notend9.dat',status='unknown',recl=512) knt = 0 ! count all tracks in last region ksrg = -1 ! flag for region change ! ! --------------------------------------------------------- ! Position for009 at first particle in the last region ! read(9,*) ! skip header records read(9,*) read(9,*) ! 100 continue read(9,*,iostat=ioc) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol ! if( ioc .gt. 0 ) then write(*,*)'Error reading FOR009.DAT' go to 900 else if( ioc .eq. -1 ) then ! EOF do i=1,knt+1 backspace(unit=9) end do go to 200 end if ! if( jsrg .ne. ksrg ) then ! new region knt = 1 ksrg = jsrg else ! same region knt = knt + 1 end if ! go to 100 ! 200 continue ! --------------------------------------------------------- ! Make list of good tracks in the last region on the file ! kk = 0 210 continue read(9,*,iostat=ioc) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol ! if( ioc .gt. 0 ) then write(*,*)'Error reading FOR009.DAT' go to 900 else if( ioc .eq. -1 ) then ! EOF go to 220 end if ! if( ipflg .ne. 0 ) go to 210 ! kk = kk + 1 igd(kk) = ievt go to 210 ! 220 continue write(*,225) kk 225 format(i6,' good tracks at end of FOR009.DAT') if( kk .gt. mxtrk ) then write(*,*) 'Too many good tracks' write(*,'(a,i8)')'Only processing ',mxtrk kk = mxtrk end if ! --------------------------------------------------------- ! Copy header lines rewind(9) ! read(9,'(a80)') line write(10,'(a80)') line read(9,'(a80)') line write(10,'(a80)') line read(9,*) ! skip line of column labels read(9,'(a80)') line ! read beginning of 1st data row as characters do i=1,80 if( line(i:i) .eq. '.' ) idec = i if( line(i:i) .eq. 'E' ) iexp = i if( line(i:i).eq.' ' .and. iexp.gt.0 ) then if( i-iexp .lt. 5 ) then f9dp = iexp - idec - 1 else f9dp = 17 end if go to 250 end if end do 250 continue if( f9dp .eq. 4 ) then write(10,254) 254 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t43,'x',t55,'y',t67,'z',t79,'Px',t91,'Py', & t103,'Pz',t115,'Bx',t127,'By',t139,'Bz',t151,'wt', & t163,'Ex',t175,'Ey',t187,'Ez',t199,'arclength',t211,'polX', & t223,'polY',t235,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e12.4)' else if( f9dp .eq. 6 ) then write(10,256) 256 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t45,'x',t59,'y',t73,'z',t87,'Px',t101,'Py', & t115,'Pz',t129,'Bx',t145,'By',t157,'Bz',t171,'wt', & t185,'Ex',t199,'Ey',t213,'Ez',t227,'arclength',t241,'polX', & t255,'polY',t269,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e14.6)' else if( f9dp .eq. 8 ) then write(10,258) 258 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t47,'x',t63,'y',t79,'z',t95,'Px',t111,'Py', & t127,'Pz',t143,'Bx',t159,'By',t175,'Bz',t191,'wt', & t207,'Ex',t223,'Ey',t239,'Ez',t255,'arclength',t271,'polX', & t287,'polY',t303,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e16.8)' else if( f9dp .eq. 10 ) then write(10,260) 260 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t49,'x',t67,'y',t85,'z',t103,'Px',t121,'Py', & t139,'Pz',t157,'Bx',t175,'By',t193,'Bz',t211,'wt', & t229,'Ex',t247,'Ey',t265,'Ez',t283,'arclength',t301,'polX', & t319,'polY',t337,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e18.10)' else if( f9dp .eq. 12 ) then write(10,262) 262 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t51,'x',t71,'y',t91,'z',t111,'Px',t131,'Py', & t151,'Pz',t171,'Bx',t191,'By',t211,'Bz',t231,'wt', & t251,'Ex',t271,'Ey',t291,'Ez',t311,'arclength',t331,'polX', & t351,'polY',t371,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e20.12)' else if( f9dp .eq. 14 ) then write(10,264) 264 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t53,'x',t75,'y',t97,'z',t119,'Px',t141,'Py', & t163,'Pz',t185,'Bx',t207,'By',t229,'Bz',t251,'wt', & t273,'Ex',t295,'Ey',t317,'Ez',t339,'arclength',t361,'polX', & t383,'polY',t405,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e22.14)' else if( f9dp .eq. 16 ) then write(10,266) 266 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t55,'x',t79,'y',t103,'z',t127,'Px',t151,'Py', & t175,'Pz',t199,'Bx',t223,'By',t247,'Bz',t271,'wt', & t295,'Ex',t319,'Ey',t343,'Ez',t367,'arclength',t391,'polX', & t415,'polY',t439,'polZ') fmt = '(i7,2i4,i5,i6,1x,1p,18e24.16)' else if( f9dp .eq. 17 ) then write(10,267) 267 format('#',t5,'evt',t9,'par',t13,'typ',t18,'flg',t24,'reg', & t31,'time',t56,'x',t81,'y',t106,'z',t131,'Px',t156,'Py', & t181,'Pz',t206,'Bx',t231,'By',t256,'Bz',t281,'wt', & t306,'Ex',t331,'Ey',t356,'Ez',t381,'arclength',t406,'polX', & t431,'polY',t456,'polZ') fmt = '(i7, 2i4, i5, i6, 1x, 1p, 18e25.16E3)' else write(10,266) fmt = '(i7,2i4,i5,i6,1x,1p,18e24.16)' end if backspace(9) ! rewind 1st data line ! ! --------------------------------------------------------- ! Position file at start of first valid region number ! Watch out for files with just a reference particle at production read(9,*,iostat=ioc) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol if( ievt.eq.0 .and. jsrg.eq.1 ) then ! ref particle read(9,*,iostat=ioc) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol if( ievt.ne.0 .and. jsrg.eq.1 ) then ! full data at production backspace(9) backspace(9) else ! just ref particle at production backspace(9) end if else ! no reference particle data backspace(9) end if ksrg = jsrg ! ! --------------------------------------------------------- ! Write new file with good tracks eliminated ! Tag and count set of bad tracks nbad = 0 ! 300 continue read(9,*,end=350) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol ! do i=1,kk if( ievt .eq. igd(i) ) then go to 300 ! skip if you find good one end if end do ! write(10,fmt) ievt,ipnum,iptyp,ipflg,jsrg, & tp,xp,pp,bfld,evtwt,efld,sarc,pol ! if( jsrg .eq. ksrg ) then nbad = nbad + 1 ibd(nbad) = ievt end if go to 300 ! ! --------------------------------------------------------- 350 continue write(*,355) nbad 355 format(i6,' bad tracks in FOR009.DAT') ! if( nbad .gt. 0 ) then write(*,*) 'Set of example bad tracks' mxt = MIN(nbad,100) do i=1,mxt,10 write(*,'(10i6)') (ibd(j),j=i,i+9) end do end if ! close(9) close(10) ! 900 continue end