*$ CREATE SOURCE.FOR *COPY SOURCE * *=== source ===========================================================* * SUBROUTINE SOURCE ( NOMORE ) INCLUDE '(DBLPRC)' INCLUDE '(DIMPAR)' INCLUDE '(IOUNIT)' * *----------------------------------------------------------------------* * * * Copyright (C) 1990-2010 by Alfredo Ferrari & Paola Sala * * All Rights Reserved. * * * * * * New source for FLUKA9x-FLUKA20xy: * * * * Created on 07 January 1990 by Alfredo Ferrari & Paola Sala * * Infn - Milan * * * * Last change on 17-Oct-10 by Alfredo Ferrari * * * * This is just an example of a possible user written source routine. * * note that the beam card still has some meaning - in the scoring the * * maximum momentum used in deciding the binning is taken from the * * beam momentum. Other beam card parameters are obsolete. * * * * Output variables: * * * * Nomore = if > 0 the run will be terminated * * * *----------------------------------------------------------------------* * INCLUDE '(BEAMCM)' INCLUDE '(FHEAVY)' INCLUDE '(FLKSTK)' INCLUDE '(IOIOCM)' INCLUDE '(LTCLCM)' INCLUDE '(PAPROP)' INCLUDE '(SOURCM)' INCLUDE '(SUMCOU)' * PARAMETER (NMAX=1000000) LOGICAL LFIRST * CHARACTER*250 LINE INTEGER NNN,IJ(NMAX) DIMENSION XXX(NMAX),YYY(NMAX),ZZZ(NMAX) DIMENSION UUU(NMAX),VVV(NMAX),WWW(NMAX) DIMENSION ERG(NMAX),WGT(NMAX) SAVE LFIRST SAVE XXX,YYY,ZZZ SAVE UUU,VVV,WWW SAVE ERG,WGT DATA LFIRST / .TRUE. / *======================================================================* * * * BASIC VERSION * * * *======================================================================* * NOMORE = 0 * +-------------------------------------------------------------------* * | First call initializations: IF ( LFIRST ) THEN * | *** The following 3 cards are mandatory *** TKESUM = ZERZER LFIRST = .FALSE. LUSSRC = .TRUE. * | *** User initialization *** LUNRD = NINT(WHASOU(1)) NNN = 0 10 CONTINUE READ( LUNRD, '(A)', ERR=10, END=20 ) LINE READ (LINE,*,ERR=10) I, X, Y, Z, U, V, W, E, WG NNN = NNN + 1 IF (NNN.GT.NMAX) CALL FLABRT('SOURCE','Increase NMAX') IJ(NNN) = I XXX(NNN) = X YYY(NNN) = Y ZZZ(NNN) = Z * Normalize direction to 1.0 UVW = SQRT(U**2 + V**2 + W**2) UUU(NNN) = U / UVW VVV(NNN) = V / UVW WWW(NNN) = W / UVW ERG(NNN) = E WGT(NNN) = WG PRINT*,"I,X,Y,Z,U,V,W,E,WG",I, X, Y, Z, U, V, W, E, WG GOTO 10 20 CONTINUE IF (NNN.EQ.0) CALL FLABRT('SOURCE','Error reading file') WRITE (LUNOUT,*) WRITE (LUNOUT,*) '*** rdsource: ',NNN,' particles loaded' WRITE (LUNOUT,*) END IF NPFLKA = NPFLKA + 1 RNDSIG = FLRNDM (RNDSIG) N = INT(NNN*RNDSIG)+1 WTFLK (NPFLKA) = WGT(N) WEIPRI = WEIPRI + WTFLK (NPFLKA) PRINT*,"...DEBUG...DEBUG...WEIPRI:",WEIPRI,WGT(N) ILOFLK (NPFLKA) = IJ(N) * From this point .... LOFLK (NPFLKA) = 1 ! Generation is 1 for source particles LOUSE (NPFLKA) = 0 ! User variables: the user can set DO 100 ISPR = 1, MKBMX1 ! different values in the STUPRF or SPAREK (1,NPFLKA) = ZERZER ! STUPRE routine, but it is better 100 CONTINUE ! not to do it here DO 200 ISPR = 1, MKBMX2 ! ISPARK (ISPR,NPFLKA) = 0 ! More user variables (integer) 200 CONTINUE NPARMA = NPARMA + 1 ! Updating the maximum particle number NUMPAR (NPFLKA) = NPARMA ! Setting the current particle number NEVENT (NPFLKA) = 0 ! Resetting the current event number DFNEAR (NPFLKA) = +ZERZER ! Resetting the distance to the ! nearest boundary * ... to this point: don't change anything AGESTK (NPFLKA) = +ZERZER ! Particle age is zero by default AKNSHR (NPFLKA) = -TWOTWO ! Resets the Kshort component of ! K0/K0bar. Usually not to be changed. IGROUP (NPFLKA) = 0 ! Group number for low energy ! neutrons: if set to 0, the program ! derives it from the kinetic energy IJBEAM = IJ(N) * +-------------------------------------------------------------------* * | (Radioactive) isotope: IF ( IJBEAM .EQ. -2 .AND. LRDBEA ) THEN IARES = IPROA IZRES = IPROZ IISRES = IPROM CALL STISBM ( IARES, IZRES, IISRES ) IJHION = IPROZ * 1000 + IPROA IJHION = IJHION * 100 + KXHEAV IONID = IJHION CALL DCDION ( IONID ) CALL SETION ( IONID ) * | * +-------------------------------------------------------------------* * | Heavy ion: ELSE IF ( IJBEAM .EQ. -2 ) THEN IJHION = IPROZ * 1000 + IPROA IJHION = IJHION * 100 + KXHEAV IONID = IJHION CALL DCDION ( IONID ) CALL SETION ( IONID ) ILOFLK (NPFLKA) = IJHION * | Flag this is prompt radiation LRADDC (NPFLKA) = .FALSE. * | Group number for "low" energy neutrons, set to 0 anyway IGROUP (NPFLKA) = 0 * | * +-------------------------------------------------------------------* * | Normal hadron: ELSE IONID = IJBEAM ILOFLK (NPFLKA) = IJBEAM * | Flag this is prompt radiation LRADDC (NPFLKA) = .FALSE. * | Group number for "low" energy neutrons, set to 0 anyway IGROUP (NPFLKA) = 0 END IF TKEFLK (NPFLKA) = ERG(N) PMOFLK (NPFLKA) = SQRT(ERG(N)*(ERG(N)+TWOTWO*AM(IONID))) PRINT*,"...DEBUG...DEBUG...",RNDSIG,"N",N,'NPFLKA',NPFLKA PRINT*,"UUU(N)",UUU(N),"VVV(N)",VVV(N),"WWW(N)",WWW(N) PRINT*,"XXX(N)",XXX(N),"YYY(N)",YYY(N),"ZZZ(N)",ZZZ(N) * Cosines (tx,ty,tz) TXFLK (NPFLKA) = UUU(N) TYFLK (NPFLKA) = VVV(N) TZFLK (NPFLKA) = WWW(N) * Particle coordinates * XFLK (NPFLKA) = XXX(N) !+ XBEAM YFLK (NPFLKA) = YYY(N) + TWOTWO*TWOTWO !+ YBEAM ZFLK (NPFLKA) = ZZZ(N) !+ ZBEAM TXPOL (NPFLKA) = -TWOTWO TYPOL (NPFLKA) = +ZERZER TZPOL (NPFLKA) = +ZERZER IF ( ILOFLK (NPFLKA) .EQ. -2 .OR. ILOFLK (NPFLKA) .GT. 100000 ) & THEN TKESUM = TKESUM + TKEFLK (NPFLKA) * WTFLK (NPFLKA) ELSE IF ( ILOFLK (NPFLKA) .NE. 0 ) THEN TKESUM = TKESUM + ( TKEFLK (NPFLKA) + AMDISC (ILOFLK(NPFLKA)) ) & * WTFLK (NPFLKA) ELSE TKESUM = TKESUM + TKEFLK (NPFLKA) * WTFLK (NPFLKA) END IF RADDLY (NPFLKA) = ZERZER * Here we ask for the region number of the hitting point. * NREG (NPFLKA) = ... * The following line makes the starting region search much more * robust if particles are starting very close to a boundary: CALL GEOCRS ( TXFLK (NPFLKA), TYFLK (NPFLKA), TZFLK (NPFLKA) ) CALL GEOREG ( XFLK (NPFLKA), YFLK (NPFLKA), ZFLK (NPFLKA), & NRGFLK(NPFLKA), IDISC ) * Do not change these cards: CALL GEOHSM ( NHSPNT (NPFLKA), 1, -11, MLATTC ) NLATTC (NPFLKA) = MLATTC CMPATH (NPFLKA) = ZERZER CALL SOEVSV RETURN *=== End of subroutine Source =========================================* END