C...Sample job to demonstrate use of LUCIFER version 2.2. C...Should be loaded with LUCIFER 2.2, JETSET 6.2 and CERN.PACKLIB COMMON /LUJETS/ N,K(2000,2),P(2000,5) COMMON /PYPARA/ IPY(30),PYPAR(20),PYVAR(10) COMMON /PYSUBS/ ISELEC,ISUBPR(4),MESHT(20),IHEL COMMON /PYCROS/ CROSS(0:4),XHTMES(0:20) COMMON /PYPROC/ IN,ICONF,JFL,MESON,MHELIC,X,SHAT,THAT,Q2 COMMON // HMEMOR(5000) DIMENSION NCONF(4),INFL(-4:4),JFLFUS(4),IHTMES(0:20), &HELMES(0:1,0:20),PTMESI(0:20) CHARACTER*10 CMES(0:20) DATA CMES/'all mesons', &'pi+','pi-','pi0','eta','eta''','K+','K-','K0','G0S', &'G0A','rho+','rho-','rho0','omega','phi','K*+','K*-','K*0', &'K*0bar','G2S'/ DATA NCONF/4*0/,NEV/0/,INFL,JFLFUS,IHTMES/34*0/,HELMES/42*0./ DATA PTMESI/21*0./ C...Selection of included processes; 0->all, 1->fusion, 2->QCD Compton, C...3->higher twist, 4->QED Compton, 5->fusion+QCD Compton (default). ISELEC=0 C...Include gluonium production (off by default). MESHT(9)=1 MESHT(10)=1 MESHT(20)=1 C...Only higher twist into a single meson, e.g. eta'. C DO 10 M=1,20 C 10 MESHT(M)=0 C MESHT(5)=1 C...Skip fragmentation, e.g. when studying only the prompt meson. C IPY(7)=0 C...Trigger pt. PTTRIG=2.0 C...Use PYTIME to monitor the amount of CPU time used. CALL PYTIME(TIME1) C...Initialize for 150 GeV photon beam on proton target; min parton pt. QTMIN=PTTRIG CALL LUCINI('FIXT','P',150.,QTMIN) CALL PYTIME(TIME2) C...Print parameter values. CALL PYPRNT CALL LULIST(5) CALL HLIMIT(5000) CALL HBOOK1(2,'dsigma/dx target parton$',40,0.,1.) CALL HBOOK1(3,'dsigma/dQ$',40,0.,10.) CALL HBOOK1(4,'dsigma/dpt charged ptcls above cut$',40,0.,10.) C...Loop over events to be generated. DO 500 IEVT=1,100000 CALL TIMEL(TIML) IF(TIML.LT.3) GOTO 600 C...Generate one event. CALL LUCIFR NEV=NEV+1 C...Count the number of different scattering subprocesses. NCONF(ICONF)=NCONF(ICONF)+1 C...List first event of each kind. IF(NCONF(ICONF).EQ.1) THEN CALL LULIST(1) WRITE(6,*) ' IN, ICONF, JFL, MESON, MHELIC = ',IN,ICONF,JFL, & MESON,MHELIC ENDIF CALL HF1(2,X,1.) CALL HF1(3,SQRT(Q2),1.) IF(IN.GE.-4.AND.IN.LE.4) INFL(IN)=INFL(IN)+1 IF(ICONF.EQ.1) THEN C...Fusion into q-qbar pair. IF(JFL.GE.1.AND.JFL.LE.4) JFLFUS(JFL)=JFLFUS(JFL)+1 ELSEIF(ICONF.EQ.3) THEN C...Higher twist process. IHTMES(0)=IHTMES(0)+1 IHTMES(MESON)=IHTMES(MESON)+1 MMH=MHELIC IF(MHELIC.LT.0.OR.MHELIC.GT.1) MMH=0 HELMES(MMH,0)=HELMES(MMH,0)+1. HELMES(MMH,MESON)=HELMES(MMH,MESON)+1. ENDIF IF(IPY(7).EQ.0.AND.ICONF.NE.3) GOTO 500 C...Loop over partons/particles. DO 100 I=1,N IF(IABS(K(I,2)).GT.100) GOTO 100 PT=SQRT(P(I,1)**2+P(I,2)**2) IF(PT.LT.PTTRIG) GOTO 100 MES=MESLU(K(I,2)) IF(MES.EQ.0) GOTO 100 PTMESI(0)=PTMESI(0)+1. PTMESI(MES)=PTMESI(MES)+1. 100 CONTINUE C...Keep only charged stable particles. CALL LUEDIT(3) DO 200 I=1,N PT=SQRT(P(I,1)**2+P(I,2)**2) C...pt spectrum of charged particles above pt-cut. IF(PT.GT.PTTRIG) CALL HF1(4,PT,1.) 200 CONTINUE 500 CONTINUE 600 CONTINUE CALL PYTIME(TIME3) CALL LUCROS IF(NEV.EQ.0) THEN WRITE(6,*) ' No event generated (NEV=0)]' STOP ENDIF WRITE(6,2100) DO 700 I=0,20 700 WRITE(6,2200) I,CMES(I),IFIX(PTMESI(I)), &100.*PTMESI(I)/MAX(1.,PTMESI(0)),PTMESI(I)*CROSS(0)/NEV*1.E+6 WRITE(6,*) WRITE(6,*) ' PTTRIG = ',PTTRIG WRITE(6,*) ' QTMIN = ',QTMIN WRITE(6,*) ' Time for initialization = ',TIME2-TIME1,' seconds.' WRITE(6,*) ' Time for event generation = ',TIME3-TIME2,' seconds.' WRITE(6,*) ' Time needed for one event = ', &(TIME3-TIME2)/MAX(1,NEV),' seconds.' WRITE(6,1200) INFL,(100.*INFL(I)/NEV,I=-4,4) WRITE(6,1300) JFLFUS,(100.*JFLFUS(I)/MAX(1,NCONF(1)),I=1,4) WRITE(6,1400) IHTMES,(100.*IHTMES(I)/MAX(1,NCONF(3)),I=0,20) WRITE(6,1500) (INT(HELMES(0,J)),J=0,20), &(100.*HELMES(0,J)/MAX(1,NCONF(3)),J=0,20) WRITE(6,1600) (INT(HELMES(1,J)),J=0,20), &(100.*HELMES(1,J)/MAX(1,NCONF(3)),J=0,20) CALL HOPERA(2,'+',2,2,40./1./NEV*CROSS(0),0.) CALL HOPERA(3,'+',3,3,40./10./NEV*CROSS(0),0.) CALL HOPERA(4,'+',4,4,40./10./NEV*CROSS(0),0.) CALL HLOGAR(0) CALL HPRINT(0) 1200 FORMAT(//,1X,'# reacting flavours:',9I6, & /,1X,'% reacting flavours:',9F6.1) 1300 FORMAT(//,1X,'# fusion flavours :',4I6, & /,1X,'% fusion flavours :',4F6.1) 1400 FORMAT(//,1X,'# meson flav.',11I6,/,20X,10I6, & /,1X,'% meson flav.',11F6.1,/,20X,10F6.1) 1500 FORMAT(//,1X,'# hel=0 meson',11I6,/,20X,10I6, & /,1X,'% hel=0 meson',11F6.1,/,20X,10F6.1) 1600 FORMAT(//,1X,'# hel=1 meson',11I6,/,20X,10I6, & /,1X,'% hel=1 meson',11F6.1,/,20X,10F6.1) 2100 FORMAT(//,5X,'Cross-sections for inclusive high-pt mesons:', & /,5X,'-------------------------------------------', &//,10X,'meson',5X,22X,'events',2X,7X,'cross-section', &/,5X,'no',3X,'type',6X,14X,'#',14X,'%',15X,'in nb',/,5X,65('-')) 2200 FORMAT(5X,I2,3X,A10,I15,F15.3,5X,1P,E15.3) END C******************************************************************** FUNCTION MESLU(I) MESLU=0 IF(I.EQ.17) THEN MESLU=1 ELSEIF(I.EQ.-17) THEN MESLU=2 ELSEIF(I.EQ.23) THEN MESLU=3 ELSEIF(I.EQ.24) THEN MESLU=4 ELSEIF(I.EQ.25) THEN MESLU=5 ELSEIF(I.EQ.18) THEN MESLU=6 ELSEIF(I.EQ.-18) THEN MESLU=7 ELSEIF(I.EQ.19) THEN MESLU=8 ELSEIF(I.EQ.91) THEN MESLU=9 ELSEIF(I.EQ.92) THEN MESLU=10 ELSEIF(I.EQ.27) THEN MESLU=11 ELSEIF(I.EQ.-27) THEN MESLU=12 ELSEIF(I.EQ.33) THEN MESLU=13 ELSEIF(I.EQ.34) THEN MESLU=14 ELSEIF(I.EQ.35) THEN MESLU=15 ELSEIF(I.EQ.28) THEN MESLU=16 ELSEIF(I.EQ.-28) THEN MESLU=17 ELSEIF(I.EQ.29) THEN MESLU=18 ELSEIF(I.EQ.-29) THEN MESLU=19 ELSEIF(I.EQ.93) THEN MESLU=20 ENDIF RETURN END