IAP GITLAB

Commit c6ba9801 authored by Colin Baus's avatar Colin Baus

dpmjet goes through decay of epos. and uses hepmcstore

git-svn-id: https://devel-ik.fzk.de/svn/mc/crmc/trunk@4126 c7a5e08c-de06-0410-9364-b41cf42a0b17
parent aba1ed62
......@@ -147,7 +147,7 @@ c Fix final particles and some event parameters
call afinal
c Fill HEP common
if(model.ne.4.and.model.ne.5.and.model.ne.12)call hepmcstore
if(model.ne.4.and.model.ne.5)call hepmcstore
c optional Statistic information (only with debug level ish=1)
call astati
......
INCLUDE_DIRECTORIES ("${PROJECT_SOURCE_DIR}/src/dpmjet/3.0-6")
FILE(GLOB files *.f)
set_source_files_properties( dpmjet3.0-6.f dpmjet-crmc.f phojet1.12-35c4.f pythia6115dpm3v1.f ../../crmc-aaa.f ../../models.F PROPERTIES COMPILE_FLAGS "-D __DPMJET__")
IF (__CRMCSTATIC__)
add_library(Dpmjet STATIC dpmjet3.0-6.f dpmjet-crmc.f phojet1.12-35c4.f pythia6115dpm3v1.f ../../crmc-aaa.f ../../models.F)
ELSE (__CRMCSTATIC__)
add_library(Dpmjet SHARED dpmjet3.0-6.f dpmjet-crmc.f phojet1.12-35c4.f pythia6115dpm3v1.f ../../crmc-aaa.f ../../models.F)
target_link_libraries(Dpmjet EposBasic)
target_link_libraries(Dpmjet CrmcBasic)
ENDIF (__CRMCSTATIC__)
INSTALL (TARGETS Dpmjet
RUNTIME DESTINATION bin
......
......@@ -3,7 +3,7 @@ c authors C. Baus, A. Feydnitch
c-----------------------------------------------------------------------
subroutine IniDPMJET
c-----------------------------------------------------------------------
c Primary initialization for PHOJET
c Primary initialization for DPMJET
c-----------------------------------------------------------------------
IMPLICIT DOUBLE PRECISION(A-H,O-Z)
include 'epos.inc'
......@@ -13,18 +13,37 @@ c-----------------------------------------------------------------------
COMMON /DTPART/ ANAME(210),AAM(210),GA(210),TAU(210),
& IICH(210),IIBAR(210),K1(210),K2(210)
INTEGER NCASES,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,IGLAU
DOUBLE PRECISION EPROJ
c change to center of mass system
c common block var iframe renamed to dframe since it is used by epos
LOGICAL LEMCCK,LHADRO,LSEADI,LEVAPO
COMMON /DTFLG1/ IFRAG(2),IRESCO,IMSHL,IRESRJ,IOULEV(6),
& LEMCCK,LHADRO(0:9),LSEADI,LEVAPO,DFRAME,ITRSPT
COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
C switch off decays
MSTJ(21)=0
DFRAME=2
C general initialization
C set some variables for decays handled by epos
egymin=10. !?
egymax=1e7 !?
irescl=0 !momentum rescale
END
c-----------------------------------------------------------------------
subroutine IniEvtDpm
c-----------------------------------------------------------------------
c Setting energy,... for each event. Usefull for e.g. conex
c-----------------------------------------------------------------------
IMPLICIT DOUBLE PRECISION(A-H,O-Z)
include 'epos.inc'
INTEGER NPMASS,NPCHAR,NTMASS,NTCHAR,IDP
DOUBLE PRECISION EPROJ
COMMON /DPMEVTINI/ EPROJ,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP
C general initialization
NCASES = -1
NPMASS = maproj
NPCHAR = laproj
......@@ -37,11 +56,10 @@ C general initialization
CALL DT_INIT(NCASES,EPROJ,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,IGLAU)
END
c-----------------------------------------------------------------------
subroutine emsdpmjet(iret)
c-----------------------------------------------------------------------
c call PHOJET to simulate interaction
c call DPMJET to simulate interaction
c-----------------------------------------------------------------------
IMPLICIT DOUBLE PRECISION(A-H,O-Z)
include 'epos.inc'
......@@ -66,27 +84,21 @@ c-----------------------------------------------------------------------
COMMON /DTPRTA/ IT,ITZ,IP,IPZ,IJPROJ,IBPROJ,IJTARG,IBTARG
INTEGER NPMASS,NPCHAR,NTMASS,NTCHAR,IDP
INTEGER ITRY, IREJ, KKMAT
INTEGER IREJ, KKMAT
DOUBLE PRECISION EPROJ
COMMON /DPMEVTINI/ EPROJ,NPMASS,NPCHAR,NTMASS,NTCHAR,IDP
integer nptlhep(NMXHKK)
C general initialization
NPMASS = maproj
NPCHAR = laproj
NTMASS = matarg
NTCHAR = latarg
IDP = 1
EPROJ = (ECMS**2-AAM(NPMASS)**2-AAM(NTMASS)**2)
& /(2.0D0*AAM(NTMASS))
call conre !! adds beam particles that are knonw from intial call of crmc
call conwr
ITRY = 0
IREJ = 0
KKMAT = -1
ITRY = ITRY+1
C call dpmjet event simulation
CALL dt_kkinc(NPMASS,NPCHAR,NTMASS,NTCHAR,IDP,EPROJ, KKMAT,IREJ)
IF(IREJ.NE.0) THEN
IF(ish.ge.2) WRITE(ifch,'(a)')
WRITE(ifch,'(a)')
&'Bad return from DPMJET. Skipping event'
GOTO 1001
ENDIF
......@@ -109,14 +121,14 @@ c jpnevt ........ number of absolute projectile neutron spectators
c jppevt ........ number of absolute projectile proton spectators
c jtnevt ........ number of absolute target neutron spectators
c jtpevt ........ number of absolute target proton spectators
c xbjevt ........ bjorken x for dis
c qsqevt ........ q**2 for dis
c xbjevt ........ bjorken x for dis
c qsqevt ........ q**2 for dis
c sigtot ........ total cross section
c nglevt ........ number of collisions acc to Glauber
c zppevt ........ average Z-parton-proj
c nglevt ........ number of collisions acc to Glauber
c zppevt ........ average Z-parton-proj
c zptevt ........ average Z-parton-targ
c ng1evt ........ number of Glauber participants with at least one IAs
c ng2evt ........ number of Glauber participants with at least two IAs
c ng1evt ........ number of Glauber participants with at least one IAs
c ng2evt ........ number of Glauber participants with at least two IAs
c ikoevt ........ number of elementary parton-parton scatterings
c typevt ........ type of event (1=Non Diff, 2=Double Diff, 3=Central Diff, 4=AB->XB, -4=AB->AX)
......@@ -133,52 +145,76 @@ c PLEASE FILL!!
phievt=0.
phi=0.
IF(NHKK.GE.NMXHEP) THEN
IF(ish.ge.2) WRITE(ifch,'(a)')
&'NMAXHEP too small for event. Skipping'
IF(NHKK.GE.NMXHKK) THEN
WRITE(ifch,'(a)')
&'NMAXHKK too small for event. Skipping'
GOTO 1001
ENDIF
NEVHEP = NEVHKK
NHEP = 0
C particle copy loop
DO I=3,NHKK
c Destroys HEPEVT structure if only beam and final (status={1,4}) particles are copied
c IF(ISTHKK(I).NE.1 .AND. ISTHKK(I).NE.4) CYCLE
NHEP = NHEP + 1
ISTHEP(I-2) = ISTHKK(I)
IDHEP(I-2) = IDHKK(I)
JMOHEP(1,I-2) = JMOHKK(1,I)-2
JMOHEP(2,I-2) = JMOHKK(2,I)-2
JDAHEP(1,I-2) = JDAHKK(1,I)-2
JDAHEP(2,I-2) = JDAHKK(2,I)-2
IF(JMOHKK(1,I).EQ.0) JMOHEP(1,I-2) = 1
IF(JMOHKK(2,I).EQ.0) JMOHEP(2,I-2) = 1
IF(JDAHKK(1,I).EQ.0) JDAHEP(1,I-2) = 0
IF(JDAHKK(2,I).EQ.0) JDAHEP(2,I-2) = 0
PHEP(1,I-2) = PHKK(1,I)
PHEP(2,I-2) = PHKK(2,I)
PHEP(3,I-2) = PHKK(3,I)
PHEP(4,I-2) = PHKK(4,I)
PHEP(5,I-2) = PHKK(5,I)
VHEP(1,I-2) = VHKK(1,I)
VHEP(2,I-2) = VHKK(2,I)
VHEP(3,I-2) = VHKK(3,I)
VHEP(4,I-2) = VHKK(4,I)
ENDDO
ISTHEP(1) = 4
ISTHEP(2) = 4
JMOHEP(1,1) = -1
JMOHEP(2,1) = -1
JMOHEP(1,2) = -1
JMOHEP(2,2) = -1
c set projectile and target as non final particles
istptl(1)=1
istptl(2)=1
do k=1,NHKK
c LIST is the code of final particle, P - its 4-momentum and mass.
ic=IDHKK(k) !! copies IDs
if(ish.ge.7)write(ifch,'(a,i5,a,i5,2a,4(e10.4,1x))')
$ ' DPMJET particle ',k,' id :',ic,' before conversion'
$ , ' momentum :',(sngl(PHKK(i,k)),i=1,5)
IF(ISTHKK(k).EQ.1)THEN !! if final particle
c final particle
nptl=nptl+1 !! add 1 particle to stack
nptlhep(k)=nptl
istptl(nptl)=ISTHKK(k)-1 !!somehow 0 is final?
else
nptlhep(k)=0
write(ifch,'(a)')
& 'Not final particle: skipped ...'
cycle
endif !! all non final particles skipped
id=idtrafo('pdg','nxs',ic) !! convert to epos ids
if(ish.ge.7)write(ifch,'(a,i5,a,i5,a)')
$ ' epos particle ',nptl,' id :',id,' after conversion'
pptl(1,nptl)=sngl(PHKK(1,k)) !P_x
pptl(2,nptl)=sngl(PHKK(2,k)) !P_y
pptl(3,nptl)=sngl(PHKK(3,k)) !P_z
pptl(4,nptl)=sngl(PHKK(4,k)) !E
pptl(5,nptl)=sngl(PHKK(5,k)) !mass
ityptl(nptl)=0
iorptl(nptl)=0
if(JMOHKK(1,k).ne.0)iorptl(nptl)=nptlhep(JMOHKK(1,k))
jorptl(nptl)=0
if(JMOHKK(2,k).ne.0)jorptl(nptl)=nptlhep(JMOHKK(2,k))
ifrptl(1,nptl)=0
ifrptl(2,nptl)=0
xorptl(1,nptl)=sngl(VHKK(1,k))*1e12
xorptl(2,nptl)=sngl(VHKK(2,k))*1e12
xorptl(3,nptl)=sngl(VHKK(3,k))*1e12
xorptl(4,nptl)=sngl(VHKK(4,k))*1e12
tivptl(1,nptl)=0.
tivptl(2,nptl)=0.
idptl(nptl)=id
write(ifch,'(a,i5,a,i5,a,i5,a,4(e10.4,1x),f6.3)')
$ ' particle from DPMJET ',nptl,' id :',idptl(nptl)
$ ,' st :',istptl(nptl)
$ , ' momentum :',(pptl(i,nptl),i=1,5)
enddo !! end of loop over particles for copying to epos block
iret=0
1000 return
1001 iret=-1
goto 1000
1001 iret=-1
goto 1000
END
......
......@@ -2905,7 +2905,7 @@ c update file names
if(model.eq.8)iversn=112 !'PHOJET '
if(model.eq.9)iversn=201125 !'FLUKA '
if(model.eq.11)iversn=300 !'QGSJETII-03'
if(model.eq.8)iversn=306 !'DPMJET '
if(model.eq.12)iversn=306 !'DPMJET '
if(model.ne.1)iverso=iversn
call IniModel(model)
endif
......
......@@ -207,6 +207,13 @@
stop'please compile with requested model'
#else
call IniEvtFlu
#endif
endif
if(model.eq.12)then
#ifndef __DPMJET__
stop'please compile with requested model'
#else
call IniEvtDpm
#endif
endif
end
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment