From: Flemming Videbaek (videbaek@sgs1.hirg.bnl.gov)
Date: Sun Dec 08 2002 - 19:13:27 EST
Hi Christian,
thanks for your comments - replys are given below.
Flemming
> Hi Flemming et al,
>
> First a question in general.
>
> Suppose you do
>
> FORTRAN/CALL gbr2c.f
> FORTRAN/CALL init_detector
> USER/CONTROL/CARDS ukin 7 input.oscar
> ...
> FORTRAN/CALL gbrfile('output.cdat')
> DO i=1,10
> USER/CONTROL/ROTANG [i]*36
> USER/CONTROL/ANALYZE gbrana 100 1 0
> ENDDO
> FORTRAN/CALL gbrend
>
> Now, suppose that the input file `input.oscar' contains exactly 100
> events. Now, it seems to me, that when the first
> `USER/CONTROL/ANALYZE' is done, then the file pointer is sitting at
> then end of `input.oscar', and when the second `USER/CONTROL/ANALYZE'
> executes the file pointer isn't rewound to the beginning - hence, the
> 2'nd and up to the 10'th `USER/CONTROL/ANALYZE' will not get any
> events from the input file `input.oscar'.
>
> Is this true? If so, is it the intended behaviour?
This was certainly the intention when I wrote the original code. The analyze
should NOT rewind/back position the input
as you point out yourself the way to get around it is to re-open the input
file as in the example. If this works for oscar I do not know, but for the
zerbra inpout files it does.
>And if so, is
> there any way to make your way around it? Something like
>
> FORTRAN/CALL gbr2c.f
> FORTRAN/CALL init_detector
> ...
> FORTRAN/CALL gbrfile('output.cdat')
> DO i=1,10
> USER/CONTROL/CARDS ukin 7 input.oscar
> USER/CONTROL/ROTANG [i]*36
> USER/CONTROL/ANALYZE gbrana 100 1 0
> ENDDO
> FORTRAN/CALL gbrend
>
> A Fortran77 question: is there anyway to force READ not to read a full
> line? For example. you want to read just one number from a line, but
> you don't want to position the file pointer at the next line. I know
> how to do it in Fortran90 (keyword ADVANCE='NO' to READ), but i
> haven't seen any thing for Fortran77. Or, is there anyway to set
> markers in the file (like C's `ftell') and then go to that mark (like
> C's `fseek')?
As far I remember, this is not possible, all assuming the mode is formatted
read.
>
> And now a BRAG issue. First a short reminder of the format of an
> OSCAR 1997a file:
>
> OSC1997a
> final_id_p_x
> <code name> <version> (<a1>, <z1>)+(<a2>, <z2>) <frame>, <beam>, <ntest>
> <event blocks>
>
> where <event block> looks like
>
> <event no> <no particles> <impact param> <reaction plane>
> <particle no> <pid> <px> <py> <pz> <e> <mass> <vx> <vy> <vz> <vt>
> ...
>
> In the subroutine `ugoscrd', there's very
> little check on the validity of the input file. For example, if
> somehow you OSCAR file is broken during write, so that a number of
> particle lines isn't written to disk, then the number of particles in
> the event header doesn't match the actual number of particle lines.
> Hence, the subroutine will read
>
> <event no> <no particles> <impact param> <reaction plane>
> <particle no> <pid> <px> <py> <pz> <e> <mass>
>
> as a particle, that is with the mappings
>
> <particle no> = <event no>
> <pid> = <no particles>
> <px> = <impact parameter>
> <py> = <reaction plane>
> <pz> = <particle no>
> <e> = <pid>
> <mass> = <px>
> <vx> = <py>
> <vy> = <pz>
> <vz> = <e>
> <vt> = <mass>
>
> and clearly the next read of an event header will read in some bogus
> numbers. Now, since Fortran will happily read an integer in the file
> into a REAL in the program, one can not even use the regular I/O
> errors to detect this case. If one could read in (or peek at) the
> <particle no> before reading in the full line, then one can easily
> check if the number corresponds to the next particle, and if it
> doesn't, then one can treat it as the event number of the next
> event. However, as I pointed out above, I don't know if it's possible
> to do that in Fortran77. So, the best thing to do, is to break out of
> `ugoscrd' if one reads an <no particles> or <impact parameter> less
> than or equal to zero.
Well one can do some checking - you can e.g. read a line into a character
variable- check the syntax of such charecter
and then read the 'formatted line
character*80 line
read(IO,fmt) line
fmt format(a)
if(line(1:n).eq.'somthing") then
skip forward or whatever
..
read(line,format) ip, px,py,px,..
along those lines. Wheter this is feasible for what you discuss I do not
jnow, but the principle of reading a line , examining it
and re-reading with a specific format is possible .
>
> Shouldn't the rotations be dealt with outside of the various `NKINE'
> cases? After all, all of the cases produces the exact same thing, a
> list of particles in the form (id, px, py, pz), and nothing else.
> Hence, the rotation stuff could be done like
>
> SUBROUTINE GUKINE
> ...
> IF (NKINE.EQ.1) THEN
> ... Read into arrays PPTL, and IDPTL, and variables IDEVT, NPTL
> ELSE IF (NKINE.EQ.1) THEN
> ... Read into arrays PPTL, and IDPTL, and variables IDEVT, NPTL
> ...
> ENDIF
>
> DO I=1,NPTL
> IF (beamdir .EQ. -1) THEN
> pptl(1,i) = -pptl(1,i)
> pptl(2,i) = -pptl(2,i)
> pptl(3,i) = -pptl(3,i)
> ENDIF
>
> pt = SQRT(pptl(1,i)**2 + pptl(2,i)**2)
> pmag = SQRT(pptl(1,i)**2 + pptl(2,i)**2 + pptl(3,i)**2)
> phi = ATAN2(pptl(1,i)**2 + pptl(2,i)**2)
> theta = ATAN2(pmag, pptl(3,i))
>
> IF (rotang .NE. 0.) phi = MOD(phi + rotang,2*pi)
>
> IF (pmag .gt. pkine(1) .and.
> $ pmag .lt. pkine(2) .and.
> $ theta .gt. pkine(3) * DEGRAD .and.
> $ theta .lt. pkine(4) * DEGRAD .and.
> $ phi .gt. pkine(5) * DEGRAD .and.
> $ phi .lt. pkine(6) * DEGRAD) THEN
> IF (rotang .NE. 0) THEN
> pptl(1,i) = pt * COS(phi)
> pptl(2,i) = pt * SIN(phi)
> ENDIF
> IF (pkine(7) .EQ. 0 .OR. pkine(7) .EQ. idptl(i) THEN
> CALL gskine(PPTL(1,i),idptl(i),nvert, 0, 0, nt)
> ENDIF
> ENDIF
> ENDDO
>
> The benefit of this approach, is that all rotations and storing is
> done in one place, and each of the different nkine cases does exactly
> one thing (with a few minor side-effects), hence making the code
> easier to read and support.
>
I think you are absolute right on this. There is an issue because not all
options do checks on a theta,p, phi range, but this can easily be coded by
setting options flag inside a given nkine if one has to make a phase space
cut or not (after the rotation).
Flemming
This archive was generated by hypermail 2.1.5 : Sun Dec 08 2002 - 19:11:46 EST