Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-05 08:21:18 

0001 C*********************************************************************
0002  
0003 C...PYSLHA
0004 C...Read/write spectrum or decay data from SLHA standard file(s).
0005 C...P. Skands
0006  
0007 C...MUPDA=1 : READ SPECTRUM ON LUN=IMSS(21)
0008 C...MUPDA=2 : LOOK FOR DECAY TABLE FOR KF=KFORIG ON LUN=IMSS(22)
0009 C...MUPDA=3 : WRITE SPECTRUM ON LUN=IMSS(23)
0010 C...(MUPDA=4 : WRITE DECAY TABLE FOR KF=KFORIG ON LUN=IMSS(24))
0011 C...MUPDA=5 : READ MASS FOR KF=KFORIG ONLY (WITH DECAY TABLE)
0012       SUBROUTINE PYSLHA(MUPDA,KFORIG,IRETRN)
0013  
0014 C...Double precision and integer declarations.
0015       IMPLICIT DOUBLE PRECISION(A-H, O-Z)
0016       IMPLICIT INTEGER(I-N)
0017       INTEGER PYK,PYCHGE,PYCOMP
0018       PARAMETER (KSUSY1=1000000,KSUSY2=2000000,KTECHN=3000000,
0019      &KEXCIT=4000000,KDIMEN=5000000)
0020 C...Commonblocks.
0021       COMMON/PYDAT1/MSTU(200),PARU(200),MSTJ(200),PARJ(200)
0022       COMMON/PYDAT2/KCHG(500,4),PMAS(500,4),PARF(2000),VCKM(4,4)
0023       COMMON/PYDAT3/MDCY(500,3),MDME(8000,2),BRAT(8000),KFDP(8000,5)
0024       COMMON/PYDAT4/CHAF(500,2)
0025       CHARACTER CHAF*16
0026       CHARACTER*40 ISAVER,VISAJE
0027       COMMON/PYINT4/MWID(500),WIDS(500,5)
0028       SAVE /PYDAT1/,/PYDAT2/,/PYDAT3/,/PYDAT4/,/PYINT4/
0029 C...SUSY blocks
0030       COMMON/PYMSSM/IMSS(0:99),RMSS(0:99)
0031       COMMON/PYSSMT/ZMIX(4,4),UMIX(2,2),VMIX(2,2),SMZ(4),SMW(2),
0032      &SFMIX(16,4),ZMIXI(4,4),UMIXI(2,2),VMIXI(2,2)
0033       COMMON/PYMSRV/RVLAM(3,3,3), RVLAMP(3,3,3), RVLAMB(3,3,3)
0034       SAVE /PYMSSM/,/PYSSMT/,/PYMSRV/
0035  
0036 C...Local arrays, character variables and data.
0037       COMMON/PYLH3P/MODSEL(200),PARMIN(100),PAREXT(200),RMSOFT(0:100),
0038      &     AU(3,3),AD(3,3),AE(3,3)
0039       COMMON/PYLH3C/CPRO(2),CVER(2)
0040       SAVE /PYLH3P/,/PYLH3C/
0041       DIMENSION MMOD(100),MSPC(100),MDEC(100)
0042 C...MMOD: flags to set for each block read in.
0043 C... 1: MODSEL     2: MINPAR     3: EXTPAR     4: SMINPUTS
0044 C...MSPC: Flags to set for each block read in.
0045 C... 1: MASS       2: NMIX       3: UMIX       4: VMIX       5: SBOTMIX
0046 C... 6: STOPMIX    7: STAUMIX    8: HMIX       9: GAUGE     10: AU
0047 C...11: AD        12: AE        13: YU        14: YD        15: YE
0048 C...16: SPINFO    17: ALPHA     18: MSOFT     19: QNUMBERS
0049       CHARACTER CPRO*12,CVER*12,CHNLIN*6
0050       CHARACTER DOC*11, CHDUM*120, CHBLCK*60
0051       CHARACTER CHINL*120,CHKF*9,CHTMP*16
0052       INTEGER VERBOS
0053       SAVE VERBOS
0054 C...Date of last Change
0055       PARAMETER (DOC='05 Mar 2007')
0056 C...MQREAD(0): Number of entries I in MQREAD
0057 C...      (I): KF code for which a QNUMBERS block has been read.
0058       DIMENSION IDC(5),KFSUSY(50),MQREAD(0:100)
0059       SAVE KFSUSY,MQREAD
0060       DATA VERBOS /1/
0061       DATA NHELLO /0/
0062       DATA KFSUSY/
0063      &1000001,1000002,1000003,1000004,1000005,1000006,
0064      &2000001,2000002,2000003,2000004,2000005,2000006,
0065      &1000011,1000012,1000013,1000014,1000015,1000016,
0066      &2000011,2000012,2000013,2000014,2000015,2000016,
0067      &1000021,1000022,1000023,1000025,1000035,1000024,
0068      &1000037,1000039,     25,     35,     36,     37,
0069      &      6,     24,     45,     46,1000045, 9*0/
0070       RMFUN(IP)=PMAS(PYCOMP(IP),1)
0071  
0072 C...Hello World
0073       IF (NHELLO.EQ.0) THEN
0074         WRITE(MSTU(11),5000) DOC
0075         NHELLO=1
0076       ENDIF
0077 
0078 C...SLHA file assumed opened by user on unit LFN, stored in IMSS(20
0079 C...+MUPDA).
0080       LFN=IMSS(20+MUPDA)
0081       IF (MUPDA.EQ.5) LFN=IMSS(21)
0082       IF (MUPDA.EQ.0) LFN=IMSS(21)
0083 C...Flag that we have not yet found whatever we were asked to find.
0084       IRETRN=1
0085  
0086 C...STOP IF LFN IS ZERO (i.e. if no LFN was given).
0087       IF (LFN.EQ.0) THEN
0088         WRITE(MSTU(11),*) '* (PYSLHA:) No valid unit given in IMSS'
0089         GOTO 9999
0090       ENDIF
0091  
0092 C...If told to read spectrum, first zero all previous information.
0093       IF (MUPDA.EQ.1) THEN
0094 C...Zero all block read flags
0095         DO 100 M=1,100
0096           MMOD(M)=0
0097           MSPC(M)=0
0098           MDEC(M)=0
0099   100   CONTINUE
0100 C...Zero all (MSSM) masses, widths, and lifetimes in PYTHIA
0101         DO 110 ISUSY=1,36
0102           KC=PYCOMP(KFSUSY(ISUSY))
0103           PMAS(KC,1)=0D0
0104           PMAS(KC,2)=0D0
0105           PMAS(KC,3)=0D0
0106           PMAS(KC,4)=0D0
0107   110   CONTINUE
0108 C...Zero all (3rd gen sfermion + gaugino/higgsino) mixing matrices.
0109         DO 130 J=1,4
0110           SFMIX(5,J) =0D0
0111           SFMIX(6,J) =0D0
0112           SFMIX(15,J)=0D0
0113           DO 120 L=1,4
0114             ZMIX(L,J) =0D0
0115             ZMIXI(L,J)=0D0
0116             IF (J.LE.2.AND.L.LE.2) THEN
0117               UMIX(L,J) =0D0
0118               UMIXI(L,J)=0D0
0119               VMIX(L,J) =0D0
0120               VMIXI(L,J)=0D0
0121             ENDIF
0122   120     CONTINUE
0123 C...Zero signed masses.
0124           SMZ(J)=0D0
0125           IF (J.LE.2) SMW(J)=0D0
0126   130   CONTINUE
0127 C...NB: RMSS array not zeroed.
0128         WRITE(MSTU(11),*)
0129      &       '* (PYSLHA:) Reading in SLHA spectrum from unit ', LFN
0130  
0131 C...If reading decays, reset PYTHIA decay counters.
0132       ELSEIF (MUPDA.EQ.2) THEN
0133         KCC=100
0134         NDC=0
0135         BRSUM=0D0
0136         DO 140 KC=1,MSTU(6)
0137           IF(KC.GT.100.AND.KCHG(KC,4).GT.100) KCC=KC
0138           NDC=MAX(NDC,MDCY(KC,2)+MDCY(KC,3)-1)
0139   140   CONTINUE
0140       ELSEIF (MUPDA.EQ.5) THEN
0141 C...Zero block read flags
0142         DO 150 M=1,100
0143           MSPC(M)=0
0144  150    CONTINUE
0145       ENDIF
0146  
0147 C............READ
0148 C...(spectrum or look for decays of KF=KFORIG or MASS of KF=KFORIG
0149       IF(MUPDA.EQ.0.OR.MUPDA.EQ.1.OR.MUPDA.EQ.2.OR.MUPDA.EQ.5) THEN
0150 C...Initialize program and version strings
0151         IF(MUPDA.EQ.1.OR.MUPDA.EQ.2) THEN
0152         CPRO(MUPDA)=' '
0153         CVER(MUPDA)=' '
0154         ENDIF
0155  
0156 C...Initialize read loop
0157         MERR=0
0158         NLINE=0
0159         CHBLCK=' '
0160 C...READ NEW LINE INTO CHINL. GOTO 300 AT END-OF-FILE.
0161   160   CHINL=' '
0162         READ(LFN,'(A120)',END=300) CHINL
0163 C...Count which line number we're at.
0164         NLINE=NLINE+1
0165         WRITE(CHNLIN,'(I6)') NLINE
0166  
0167 C...Skip comment and empty lines without processing.
0168         IF (CHINL(1:1).EQ.'#'.OR.CHINL.EQ.' ') GOTO 160
0169  
0170 C...We assume all upper case below. Rewrite CHINL to all upper case.
0171         INL=0
0172         IGOOD=0
0173   170   INL=INL+1
0174         IF (CHINL(INL:INL).NE.'#') THEN
0175           DO 180 ICH=97,122
0176             IF (CHAR(ICH).EQ.CHINL(INL:INL)) CHINL(INL:INL)=CHAR(ICH-32)
0177   180     CONTINUE
0178 C...Extra safety. Chek for sensible input on line
0179           IF (IGOOD.EQ.0) THEN
0180             DO 190 ICH=48,90
0181               IF (CHAR(ICH).EQ.CHINL(INL:INL)) IGOOD=1
0182   190       CONTINUE
0183           ENDIF
0184           IF (INL.LT.120) GOTO 170
0185         ENDIF
0186         IF (IGOOD.EQ.0) GOTO 160
0187  
0188 C...Check for BLOCK begin statement (spectrum).
0189         IF (CHINL(1:1).EQ.'B') THEN
0190           MERR=0
0191           READ(CHINL,'(A6,A)',ERR=460) CHDUM,CHBLCK
0192 C...Check if another of this type of block was already read.
0193 C...(logarithmic interpolation not yet implemented, so duplicates always
0194 C...give errors)
0195           IF (CHBLCK(1:6).EQ.'MODSEL'.AND.MMOD(1).NE.0) MERR=7
0196           IF (CHBLCK(1:6).EQ.'MINPAR'.AND.MMOD(2).NE.0) MERR=7
0197           IF (CHBLCK(1:6).EQ.'EXTPAR'.AND.MMOD(3).NE.0) MERR=7
0198           IF (CHBLCK(1:8).EQ.'SMINPUTS'.AND.MMOD(4).NE.0) MERR=7
0199           IF (CHBLCK(1:4).EQ.'MASS'.AND.MSPC(1).NE.0) MERR=7
0200           IF (CHBLCK(1:4).EQ.'NMIX'.AND.MSPC(2).NE.0) MERR=7
0201           IF (CHBLCK(1:4).EQ.'UMIX'.AND.MSPC(3).NE.0) MERR=7
0202           IF (CHBLCK(1:4).EQ.'VMIX'.AND.MSPC(4).NE.0) MERR=7
0203           IF (CHBLCK(1:7).EQ.'SBOTMIX'.AND.MSPC(5).NE.0) MERR=7
0204           IF (CHBLCK(1:7).EQ.'STOPMIX'.AND.MSPC(6).NE.0) MERR=7
0205           IF (CHBLCK(1:7).EQ.'STAUMIX'.AND.MSPC(7).NE.0) MERR=7
0206           IF (CHBLCK(1:4).EQ.'HMIX'.AND.MSPC(8).NE.0) MERR=7
0207           IF (CHBLCK(1:5).EQ.'ALPHA'.AND.MSPC(17).NE.0) MERR=7
0208           IF (CHBLCK(1:5).EQ.'AU'.AND.MSPC(10).NE.0) MERR=7
0209           IF (CHBLCK(1:5).EQ.'AD'.AND.MSPC(11).NE.0) MERR=7
0210           IF (CHBLCK(1:5).EQ.'AE'.AND.MSPC(12).NE.0) MERR=7
0211           IF (CHBLCK(1:5).EQ.'MSOFT'.AND.MSPC(18).NE.0) MERR=7
0212 C...Check for new particles
0213           IF (CHBLCK(1:8).EQ.'QNUMBERS'.OR.CHBLCK(1:8).EQ.'PARTICLE')
0214      &        THEN
0215             MSPC(19)=MSPC(19)+1
0216 C...Read PDG code
0217             READ(CHBLCK(9:60),*) KFQ
0218 
0219             DO 121 MQ=1,MQREAD(0)
0220               IF (MQREAD(MQ).EQ.KFQ) THEN
0221                 MERR=17
0222                 GOTO 290
0223               ENDIF
0224  121        CONTINUE
0225             WRITE(MSTU(11),'(A,I9,A,F12.3)')
0226      &           ' * (PYSLHA:) Reading in '//CHBLCK(1:8)//
0227      &           ' for KF =',KFQ
0228             MQREAD(0)=MQREAD(0)+1
0229             MQREAD(MQREAD(0))=KFQ
0230             MSPC(19)=MSPC(19)+1
0231             KCQ=PYCOMP(KFQ)
0232             IF (KCQ.EQ.0) THEN
0233               DO 123 KCT=100,MSTU(6)
0234                 IF(KCHG(KCT,4).GT.100) KCQ=KCT
0235  123          CONTINUE
0236               KCQ=KCQ+1
0237               KCC=KCQ
0238               KCHG(KCQ,4)=KFQ              
0239 C...First write PDG code as name
0240               WRITE(CHTMP,*) KFQ
0241 C...Then look for real name
0242               ICMT=9
0243  90           ICMT=ICMT+1
0244               IF (CHBLCK(ICMT:ICMT).NE.'#'.AND.ICMT.LT.59) GOTO 90
0245               IF (ICMT.LT.59) THEN
0246                 READ(CHBLCK(ICMT+1:60),'(A)',ERR=95) CHDUM
0247                 IF (CHDUM.NE.' ') CHTMP=CHDUM
0248               ENDIF
0249  95           IF (CHTMP(1:1).EQ.' ') THEN
0250                 READ(CHTMP,'(1x,A)') CHAF(KCQ,1)
0251               ELSE
0252                 READ(CHTMP,'(A)') CHAF(KCQ,1)
0253               ENDIF
0254               MSTU(20)=0
0255 C...Set stable for now
0256               PMAS(KCQ,2)=1D-6
0257               MWID(KCQ)=0
0258               MDCY(KCQ,1)=0
0259               MDCY(KCQ,2)=0
0260               MDCY(KCQ,3)=0
0261             ELSE
0262               WRITE(MSTU(11),*)
0263      &           '* (PYSLHA:) KF =',KFQ,' already exists: ',
0264      &             CHAF(KCQ,1), '. Entry ignored.'
0265               MERR=7
0266             ENDIF
0267           ENDIF
0268 C...Finalize this line and read next.
0269           GOTO 290
0270 C...Check for DECAY begin statement (decays).
0271         ELSEIF (CHINL(1:1).EQ.'D') THEN
0272           MERR=0
0273           BRSUM=0D0
0274           CHBLCK='DECAY'
0275 C...Read KF code and WIDTH
0276           MPSIGN=1
0277           READ(CHINL(7:INL),*,ERR=470) KF, WIDTH
0278           IF (KF.LE.0) THEN
0279             KF=-KF
0280             MPSIGN=-1
0281           ENDIF
0282 C...If this is not the KF we're looking for...
0283           IF (KF.NE.KFORIG.OR.MUPDA.NE.2) THEN
0284 C...Set block skip flag and read next line.
0285             MERR=16
0286             GOTO 290
0287           ENDIF
0288  
0289 C...Determine PYTHIA KC code of particle
0290           KCREP=0
0291           IF(KF.LE.100) THEN
0292             KCREP=KF
0293           ELSE
0294             DO 200 KCR=101,KCC
0295               IF(KCHG(KCR,4).EQ.KF) KCREP=KCR
0296   200       CONTINUE
0297           ENDIF
0298           KC=KCREP
0299           IF (KCREP.NE.0) THEN
0300 C...Particle is already known. Don't do anything yet.
0301           ELSE
0302 C...  Add new particle. Actually, this should not happen.
0303 C...  New particles should be added already when reading the spectrum
0304 C...  information, so go under previously stable category.
0305             KCC=KCC+1
0306             KC=KCC
0307           ENDIF
0308  
0309           IF (WIDTH.LE.0D0) THEN
0310 C...Stable (i.e. LSP)
0311             WRITE(MSTU(11),*)
0312      &           '* (PYSLHA:) Reading in SLHA stable particle: ',
0313      &           CHAF(KCREP,1)
0314             IF (WIDTH.LT.0D0) THEN
0315               CALL PYERRM(19,'(PYSLHA:) Negative width forced to'//
0316      &             ' zero !')
0317               WIDTH=0D0
0318             ENDIF
0319             PMAS(KC,2)=1D-6
0320             MWID(KC)=0
0321             MDCY(KC,1)=0
0322 C...Ignore any decay lines that may be present for this KF
0323             MERR=16
0324             MDCY(KC,2)=0
0325             MDCY(KC,3)=0
0326 C...Return ok
0327             IRETRN=0
0328           ENDIF
0329 C...Finalize and start reading in decay modes.
0330           GOTO 290
0331         ELSEIF (MOD(MERR,10).GE.6) THEN
0332 C...If ignore block flag set, skip directly to next line.
0333           GOTO 160
0334         ENDIF
0335  
0336 C...READ SPECTRUM
0337         IF (MUPDA.EQ.0.AND.MERR.EQ.0) THEN
0338           IF (CHBLCK(1:8).EQ.'QNUMBERS'.OR.CHBLCK(1:8).EQ.'PARTICLE') 
0339      &        THEN
0340             READ(CHINL,*) INDX, IVAL
0341             IF (INDX.EQ.1) KCHG(KCQ,1)=IVAL
0342             IF (INDX.EQ.3) KCHG(KCQ,2)=0
0343             IF (INDX.EQ.3.AND.IVAL.EQ.3) KCHG(KCQ,2)=1
0344             IF (INDX.EQ.3.AND.IVAL.EQ.-3) KCHG(KCQ,2)=-1
0345             IF (INDX.EQ.3.AND.IVAL.EQ.8) KCHG(KCQ,2)=2
0346             IF (INDX.EQ.4) THEN
0347               KCHG(KCQ,3)=IVAL
0348               IF (IVAL.EQ.1) THEN 
0349                 CHTMP=CHAF(KCQ,1)
0350                 IF (CHTMP.EQ.' ') THEN
0351                   WRITE(CHAF(KCQ,1),*) KCHG(KCQ,4)
0352                   WRITE(CHAF(KCQ,2),*) -KCHG(KCQ,4)
0353                 ELSE
0354                   ILAST=17
0355  116              ILAST=ILAST-1
0356                   IF (CHTMP(ILAST:ILAST).EQ.' ') GOTO 116
0357                   IF (CHTMP(ILAST:ILAST).EQ.'+') THEN
0358                     CHTMP(ILAST:ILAST)='-'
0359                   ELSE
0360                     CHTMP(ILAST+1:MIN(16,ILAST+4))='bar'
0361                   ENDIF
0362                   CHAF(KCQ,2)=CHTMP
0363                 ENDIF
0364               ENDIF
0365             ENDIF
0366           ELSE
0367             MERR=8
0368           ENDIF
0369         ELSEIF ((MUPDA.EQ.1.OR.MUPDA.EQ.5).AND.MERR.EQ.0) THEN
0370 C...MASS: Mass spectrum
0371           IF (CHBLCK(1:4).EQ.'MASS') THEN
0372             READ(CHINL,*) KF, VAL
0373             MERR=1
0374             KC=0
0375             IF (MUPDA.EQ.1.OR.KF.EQ.KFORIG) THEN
0376 C...Read in masses for anything
0377               MERR=0
0378               KC=PYCOMP(KF)
0379               IF (KC.NE.0) THEN
0380                 MSPC(1)=MSPC(1)+1
0381                 PMAS(KC,1) = ABS(VAL)
0382                 IF (MUPDA.EQ.5) THEN
0383                   WRITE(MSTU(11),'(A,I9,A,F12.3)')
0384      &                 ' * (PYSLHA:) Reading in MASS entry for KF =',
0385      &                 KF, ', pole mass =', VAL
0386                   IRETRN=0
0387                 ENDIF
0388 C...  Signed masses
0389                 IF (KF.EQ.1000021.AND.MSPC(18).EQ.0) RMSS(3)=VAL
0390                 IF (KF.EQ.1000022) SMZ(1)=VAL
0391                 IF (KF.EQ.1000023) SMZ(2)=VAL
0392                 IF (KF.EQ.1000025) SMZ(3)=VAL
0393                 IF (KF.EQ.1000035) SMZ(4)=VAL
0394                 IF (KF.EQ.1000024) SMW(1)=VAL
0395                 IF (KF.EQ.1000037) SMW(2)=VAL
0396               ENDIF
0397             ELSEIF (MUPDA.EQ.5) THEN
0398               MERR=0
0399             ENDIF
0400           ELSEIF (MUPDA.EQ.5) THEN
0401 C...Only read MASS if MUPDA = 5. Skip any other blocks.
0402             MERR=8
0403           ELSEIF (CHBLCK(1:8).EQ.'QNUMBERS'.OR.
0404      &          CHBLCK(1:8).EQ.'PARTICLE') THEN
0405 C...Don't print a warning for QNUMBERS when reading spectrum
0406             MERR=8
0407 C...  MODSEL: Model selection and global switches
0408           ELSEIF (CHBLCK(1:6).EQ.'MODSEL') THEN
0409             READ(CHINL,*) INDX, IVAL
0410             IF (INDX.LE.200.AND.INDX.GT.0) THEN
0411               MODSEL(INDX)=IVAL
0412               MMOD(1)=MMOD(1)+1
0413               IF (INDX.EQ.3.AND.IVAL.EQ.1) THEN
0414 C...  Switch on NMSSM
0415                 WRITE(MSTU(11),*) '* (PYSLHA:) switching on NMSSM'
0416                 IMSS(13)=MAX(1,IMSS(13))
0417 C...  Add NMSSM states if not already done
0418  
0419                 KFN=25
0420                 KCN=KFN
0421                 CHAF(KCN,1)='H_10'
0422                 CHAF(KCN,2)=' '
0423  
0424                 KFN=35
0425                 KCN=KFN
0426                 CHAF(KCN,1)='H_20'
0427                 CHAF(KCN,2)=' '
0428  
0429                 KFN=45
0430                 KCN=KFN
0431                 CHAF(KCN,1)='H_30'
0432                 CHAF(KCN,2)=' '
0433  
0434                 KFN=36
0435                 KCN=KFN
0436                 CHAF(KCN,1)='A_10'
0437                 CHAF(KCN,2)=' '
0438  
0439                 KFN=46
0440                 KCN=KFN
0441                 CHAF(KCN,1)='A_20'
0442                 CHAF(KCN,2)=' '
0443  
0444                 KFN=1000045
0445                 KCN=PYCOMP(KFN)
0446                 IF (KCN.EQ.0) THEN
0447                   DO 234 KCT=100,MSTU(6)
0448                     IF(KCHG(KCT,4).GT.100) KCN=KCT
0449  234              CONTINUE
0450                   KCN=KCN+1
0451                   KCHG(KCN,4)=KFN
0452                   MSTU(20)=0
0453                 ENDIF
0454 C...  Set stable for now
0455                 PMAS(KCN,2)=1D-6
0456                 MWID(KCN)=0
0457                 MDCY(KCN,1)=0
0458                 MDCY(KCN,2)=0
0459                 MDCY(KCN,3)=0
0460                 CHAF(KCN,1)='~chi_50'
0461                 CHAF(KCN,2)=' '
0462               ENDIF
0463             ELSE
0464               MERR=1
0465             ENDIF
0466 C...MINPAR: Minimal model parameters
0467           ELSEIF (CHBLCK(1:6).EQ.'MINPAR') THEN
0468             IF (MODSEL(1).NE.0) THEN
0469               READ(CHINL,*) INDX, VAL
0470               IF (INDX.LE.100.AND.INDX.GT.0) THEN
0471                 PARMIN(INDX)=VAL
0472                 MMOD(2)=MMOD(2)+1
0473               ELSE
0474                 MERR=1
0475               ENDIF
0476             ELSEIF (MMOD(3).NE.0) THEN
0477               WRITE(MSTU(11),*)
0478      &             '* (PYSLHA:) MINPAR after EXTPAR !'
0479               MERR=1
0480             ELSE
0481               WRITE(MSTU(11),*)
0482      &             '* (PYSLHA:) Reading MINPAR, but no MODSEL !' 
0483               MERR=1
0484             ENDIF
0485 C...tan(beta)
0486             IF (INDX.EQ.3) RMSS(5)=VAL
0487 C...EXTPAR: non-minimal model parameters.
0488           ELSEIF (CHBLCK(1:6).EQ.'EXTPAR') THEN
0489             IF (MMOD(1).NE.0) THEN
0490               READ(CHINL,*) INDX, VAL
0491               IF (INDX.LE.200.AND.INDX.GT.0) THEN
0492                 PAREXT(INDX)=VAL
0493                 MMOD(3)=MMOD(3)+1
0494               ELSE
0495                 MERR=1
0496               ENDIF
0497             ELSE
0498               WRITE(MSTU(11),*)
0499      &             '* (PYSLHA:) Reading EXTPAR, but no MODSEL !'
0500               MERR=1
0501             ENDIF
0502 C...tan(beta)
0503             IF (INDX.EQ.25) RMSS(5)=VAL
0504           ELSEIF (CHBLCK(1:8).EQ.'SMINPUTS') THEN
0505             READ(CHINL,*) INDX, VAL
0506             IF (INDX.LE.3.OR.INDX.EQ.5.OR.INDX.GE.7) THEN
0507               MERR=1
0508             ELSEIF (INDX.EQ.4) THEN
0509               PMAS(PYCOMP(23),1)=VAL
0510             ELSEIF (INDX.EQ.6) THEN
0511               PMAS(PYCOMP(6),1)=VAL
0512             ENDIF
0513           ELSEIF (CHBLCK(1:4).EQ.'NMIX'.OR.CHBLCK(1:4).EQ.'VMIX'.OR
0514      $           .CHBLCK(1:4).EQ.'UMIX'.OR.CHBLCK(1:7).EQ.'STOPMIX'.OR
0515      $           .CHBLCK(1:7).EQ.'SBOTMIX'.OR.CHBLCK(1:7).EQ.'STAUMIX')
0516      $           THEN
0517 C...NMIX,UMIX,VMIX,STOPMIX,SBOTMIX, and STAUMIX. Mixing.
0518             IM=0
0519             IF (CHBLCK(5:6).EQ.'IM') IM=1
0520   250       READ(CHINL,*) INDX1, INDX2, VAL
0521             IF (CHBLCK(1:1).EQ.'N'.AND.INDX1.LE.4.AND.INDX2.LE.4) THEN
0522               IF (IM.EQ.0) ZMIX(INDX1,INDX2) = VAL
0523               IF (IM.EQ.1) ZMIXI(INDX1,INDX2)= VAL
0524               MSPC(2)=MSPC(2)+1
0525             ELSEIF (CHBLCK(1:1).EQ.'U') THEN
0526               IF (IM.EQ.0) UMIX(INDX1,INDX2) = VAL
0527               IF (IM.EQ.1) UMIXI(INDX1,INDX2)= VAL
0528               MSPC(3)=MSPC(3)+1
0529             ELSEIF (CHBLCK(1:1).EQ.'V') THEN
0530               IF (IM.EQ.0) VMIX(INDX1,INDX2) = VAL
0531               IF (IM.EQ.1) VMIXI(INDX1,INDX2)= VAL
0532               MSPC(4)=MSPC(4)+1
0533             ELSEIF (CHBLCK(1:4).EQ.'STOP'.OR.CHBLCK(1:4).EQ.'SBOT'.OR
0534      $             .CHBLCK(1:4).EQ.'STAU') THEN
0535               IF (CHBLCK(1:4).EQ.'STOP') THEN
0536                 KFSM=6
0537                 ISPC=6
0538               ELSEIF (CHBLCK(1:4).EQ.'SBOT') THEN
0539                 KFSM=5
0540                 ISPC=5
0541               ELSEIF (CHBLCK(1:4).EQ.'STAU') THEN
0542                 KFSM=15
0543                 ISPC=7
0544               ENDIF
0545 C...Set SFMIX element
0546               SFMIX(KFSM,2*(INDX1-1)+INDX2)=VAL
0547               MSPC(ISPC)=MSPC(ISPC)+1
0548             ENDIF
0549 C...Running parameters
0550           ELSEIF (CHBLCK(1:4).EQ.'HMIX') THEN
0551             READ(CHBLCK(8:25),*,ERR=510) Q
0552             READ(CHINL,*) INDX, VAL
0553             MSPC(8)=MSPC(8)+1
0554             IF (INDX.EQ.1) THEN
0555               RMSS(4) = VAL
0556             ELSE
0557               MERR=1
0558               MSPC(8)=MSPC(8)-1
0559             ENDIF
0560           ELSEIF (CHBLCK(1:5).EQ.'ALPHA') THEN
0561             READ(CHINL,*,ERR=520) VAL
0562             RMSS(18)= VAL
0563             MSPC(17)=MSPC(17)+1
0564 C...Higgs parameters set manually or with FeynHiggs.
0565             IMSS(4)=MAX(2,IMSS(4))
0566           ELSEIF (CHBLCK(1:2).EQ.'AU'.OR.CHBLCK(1:2).EQ.'AD'.OR
0567      &           .CHBLCK(1:2).EQ.'AE') THEN
0568             READ(CHBLCK(9:26),*,ERR=510) Q
0569             READ(CHINL,*) INDX1, INDX2, VAL
0570             IF (CHBLCK(2:2).EQ.'U') THEN
0571               AU(INDX1,INDX2)=VAL
0572               IF (INDX1.EQ.3.AND.INDX2.EQ.3) RMSS(16)=VAL
0573               MSPC(11)=MSPC(11)+1
0574             ELSEIF (CHBLCK(2:2).EQ.'D') THEN
0575               AD(INDX1,INDX2)=VAL
0576               IF (INDX1.EQ.3.AND.INDX2.EQ.3) RMSS(15)=VAL
0577               MSPC(10)=MSPC(10)+1
0578             ELSEIF (CHBLCK(2:2).EQ.'E') THEN
0579               AE(INDX1,INDX2)=VAL
0580               IF (INDX1.EQ.3.AND.INDX2.EQ.3) RMSS(17)=VAL
0581               MSPC(12)=MSPC(12)+1
0582             ELSE
0583               MERR=1
0584             ENDIF
0585           ELSEIF (CHBLCK(1:5).EQ.'MSOFT') THEN
0586             IF (MSPC(18).EQ.0) THEN
0587               READ(CHBLCK(9:25),*,ERR=510) Q
0588               RMSOFT(0)=Q
0589             ENDIF
0590             READ(CHINL,*) INDX, VAL
0591             RMSOFT(INDX)=VAL
0592             MSPC(18)=MSPC(18)+1
0593           ELSEIF (CHBLCK(1:5).EQ.'GAUGE') THEN
0594             MERR=8
0595           ELSEIF (CHBLCK(1:2).EQ.'YU'.OR.CHBLCK(1:2).EQ.'YD'.OR
0596      &           .CHBLCK(1:2).EQ.'YE') THEN
0597             MERR=8
0598           ELSEIF (CHBLCK(1:6).EQ.'SPINFO') THEN
0599             READ(CHINL(1:6),*) INDX
0600             IT=0
0601             MIRD=0
0602   260       IT=IT+1
0603             IF (CHINL(IT:IT).EQ.' ') GOTO 260
0604 C...Don't read index
0605             IF (CHINL(IT:IT).EQ.CHAR(INDX+48).AND.MIRD.EQ.0) THEN
0606               MIRD=1
0607               GOTO 260
0608             ENDIF
0609             IF (INDX.EQ.1) CPRO(1)=CHINL(IT:IT+12)
0610             IF (INDX.EQ.2) CVER(1)=CHINL(IT:IT+12)
0611           ELSE
0612 C...  Set unrecognized block flag.
0613             MERR=6
0614           ENDIF
0615  
0616 C...DECAY TABLES
0617 C...Read in decay information
0618         ELSEIF (MUPDA.EQ.2.AND.MERR.EQ.0) THEN
0619 C...Read new decay chanel
0620           IF(CHINL(1:1).EQ.' '.AND.CHBLCK(1:5).EQ.'DECAY') THEN
0621             NDC=NDC+1
0622 C...Read in branching ratio and number of daughters for this mode.
0623             READ(CHINL(4:50),*,ERR=480) BRAT(NDC)
0624             READ(CHINL(4:50),*,ERR=490) DUM, NDA
0625             IF (NDA.LE.5) THEN
0626               IF(NDC.GT.MSTU(7)) CALL PYERRM(27,
0627      &             '(PYSLHA:) Decay data arrays full by KF ='
0628      $             //CHAF(KC,1))
0629 C...If first decay chanel, set decays start point in decay table
0630               IF(BRSUM.LE.0D0.AND.BRAT(NDC).NE.0D0) THEN 
0631                 WRITE(MSTU(11),*)
0632      &              '* (PYSLHA:) Reading in SLHA decay table for ',
0633      &              CHAF(KCREP,1)
0634 C...Set particle parameters (mass set when reading BLOCK MASS above)
0635                 PMAS(KC,2)=WIDTH
0636                 IF (KF.EQ.25.OR.KF.EQ.35.OR.KF.EQ.36) THEN
0637                   WRITE(MSTU(11),*)
0638      &                '*  Note: the Pythia gg->h/H/A cross section'//
0639      &                ' is proportional to the h/H/A->gg width'
0640                 ENDIF
0641                 PMAS(KC,3)=0D0
0642                 PMAS(KC,4)=PARU(3)*1D-12/WIDTH
0643                 MWID(KC)=2
0644                 MDCY(KC,1)=1
0645                 MDCY(KC,2)=NDC
0646                 MDCY(KC,3)=0
0647 C...Return ok
0648                 IRETRN=0
0649               ENDIF
0650 C...  Count up number of decay modes for this particle
0651               MDCY(KC,3)=MDCY(KC,3)+1
0652 C...  Read in decay daughters.
0653               READ(CHINL(4:120),*,ERR=500) DUM,IDM, (IDC(IDA),IDA=1,NDA)
0654 C...  Flip sign if reading antiparticle decays (if antipartner exists)
0655               DO 270 IDA=1,NDA
0656                 IF (KCHG(PYCOMP(IDC(IDA)),3).NE.0)
0657      &               IDC(IDA)=MPSIGN*IDC(IDA)
0658   270         CONTINUE
0659 C...Switch on decay channel, with products ordered in decreasing ABS(KF)
0660               MDME(NDC,1)=1
0661               IF (BRAT(NDC).LE.0D0) MDME(NDC,1)=0
0662               BRSUM=BRSUM+ABS(BRAT(NDC))
0663               BRAT(NDC)=ABS(BRAT(NDC))
0664  274          IFLIP=0
0665               DO 277 IDA=1,NDA-1
0666                 IF (IABS(IDC(IDA+1)).GT.IABS(IDC(IDA))) THEN
0667                   ITMP=IDC(IDA)
0668                   IDC(IDA)=IDC(IDA+1)
0669                   IDC(IDA+1)=ITMP
0670                   IFLIP=IFLIP+1
0671                 ENDIF
0672  277          CONTINUE
0673               IF (IFLIP.GT.0) GOTO 274
0674 C              WRITE(MSTU(11),7510) BRAT(NDC), NDA, (IDC(IDA),IDA=1,NDA)
0675 C...Treat as ordinary decay, no fancy stuff.
0676               MDME(NDC,2)=0
0677               DO 280 IDA=1,5
0678                 IF (IDA.LE.NDA) THEN
0679                   KFDP(NDC,IDA)=IDC(IDA)
0680                 ELSE
0681                   KFDP(NDC,IDA)=0
0682                 ENDIF
0683   280         CONTINUE
0684             ELSE
0685               CALL PYERRM(7,'(PYSLHA:) Too many daughters on line '//
0686      &             CHNLIN)
0687               MERR=11
0688               NDC=NDC-1
0689             ENDIF
0690           ELSEIF(CHINL(1:1).EQ.'+') THEN
0691             MERR=11
0692           ELSEIF(CHBLCK(1:6).EQ.'DCINFO') THEN
0693             MERR=16
0694           ELSE
0695             MERR=16
0696           ENDIF
0697         ENDIF
0698 C...  Error check.
0699   290   IF (MOD(MERR,10).EQ.1.AND.(MUPDA.EQ.1.OR.MUPDA.EQ.2)) THEN
0700           WRITE(MSTU(11),*) '* (PYSLHA:) Ignoring line '//CHNLIN//': '
0701      &         //CHINL(1:40)
0702           MERR=0
0703         ELSEIF (MERR.EQ.6.AND.MUPDA.EQ.1) THEN
0704           WRITE(MSTU(11),*) '* (PYSLHA:) Ignoring BLOCK '//
0705      &         CHBLCK(1:INL)//'... on line'//CHNLIN
0706         ELSEIF (MERR.EQ.8.AND.MUPDA.EQ.1) THEN
0707           WRITE(MSTU(11),*) '* (PYSLHA:) PYTHIA will not use BLOCK '
0708      &         //CHBLCK(1:INL)//'... on line'//CHNLIN
0709         ELSEIF (MERR.EQ.16.AND.MUPDA.EQ.2.AND.IMSS(21).EQ.0.AND.
0710      &         CHBLCK(1:1).NE.'D'.AND.VERBOS.EQ.1) THEN
0711           WRITE(MSTU(11),*) '* (PYSLHA:) Ignoring BLOCK '//CHBLCK(1:INL)
0712      &         //'... on line'//CHNLIN
0713         ELSEIF (MERR.EQ.7.AND.MUPDA.EQ.1) THEN
0714           WRITE(MSTU(11),*) '* (PYSLHA:) Ignoring extra BLOCK '/
0715      &         /CHBLCK(1:INL)//'... on line'//CHNLIN
0716         ELSEIF (MERR.EQ.2.AND.MUPDA.EQ.1) THEN
0717           WRITE (CHTMP,*) KF
0718           WRITE(MSTU(11),*)
0719      &         '* (PYSLHA:) Ignoring extra MASS entry for KF='//
0720      &         CHTMP(1:9)//' on line'//CHNLIN
0721         ENDIF
0722 C...  End of loop
0723         GOTO 160
0724   300   CONTINUE
0725 C...Set flag that KC codes have been rearranged.
0726         MSTU(20)=0
0727         VERBOS=0
0728  
0729 C...Perform possible tests that new information is consistent.
0730         IF (MUPDA.EQ.1) THEN
0731           MSTU23=MSTU(23)
0732           MSTU27=MSTU(27)
0733 C...Check Z and top masses
0734           IF (ABS(PMAS(PYCOMP(23),1)-91.2D0).GT.1D0) THEN
0735             WRITE(CHTMP,*) PMAS(PYCOMP(23),1)
0736             CALL PYERRM(19,'(PYSLHA:) note Z boson mass, M ='//CHTMP)
0737           ENDIF
0738           IF (ABS(PMAS(PYCOMP(6),1)-175D0).GT.25D0) THEN
0739             WRITE(CHTMP,*) PMAS(PYCOMP(6),1)
0740             CALL PYERRM(19,'(PYSLHA:) note top quark mass, M ='
0741      &           //CHTMP//'GeV')
0742           ENDIF
0743 C...Check masses
0744           DO 310 ISUSY=1,37
0745             KF=KFSUSY(ISUSY)
0746 C...Don't complain about right-handed neutrinos
0747             IF (KF.EQ.KSUSY2+12.OR.KF.EQ.KSUSY2+14.OR.KF.EQ.KSUSY2
0748      &           +16) GOTO 310
0749 C...Only check gravitino in GMSB scenarios
0750             IF (MODSEL(1).NE.2.AND.KF.EQ.KSUSY1+39) GOTO 310
0751             KC=PYCOMP(KF)
0752             IF (PMAS(KC,1).EQ.0D0) THEN
0753               WRITE(CHTMP,*) KF
0754               CALL PYERRM(9
0755      &             ,'(PYSLHA:) No mass information found for KF = '
0756      &             //CHTMP)
0757             ENDIF
0758   310     CONTINUE
0759 C...Check mixing matrices (MSSM only)
0760           IF (IMSS(13).EQ.0) THEN
0761             IF (MSPC(2).NE.16.AND.MSPC(2).NE.32) CALL PYERRM(9
0762      &           ,'(PYSLHA:) Inconsistent # of elements in NMIX')
0763             IF (MSPC(3).NE.4.AND.MSPC(3).NE.8) CALL PYERRM(9
0764      &           ,'(PYSLHA:) Inconsistent # of elements in UMIX')
0765             IF (MSPC(4).NE.4.AND.MSPC(4).NE.8) CALL PYERRM(9
0766      &           ,'(PYSLHA:) Inconsistent # of elements in VMIX')
0767             IF (MSPC(5).NE.4) CALL PYERRM(9
0768      &           ,'(PYSLHA:) Inconsistent # of elements in SBOTMIX')
0769             IF (MSPC(6).NE.4) CALL PYERRM(9
0770      &           ,'(PYSLHA:) Inconsistent # of elements in STOPMIX')
0771             IF (MSPC(7).NE.4) CALL PYERRM(9
0772      &           ,'(PYSLHA:) Inconsistent # of elements in STAUMIX')
0773             IF (MSPC(8).LT.1) CALL PYERRM(9
0774      &           ,'(PYSLHA:) Too few elements in HMIX')
0775             IF (MSPC(10).EQ.0) CALL PYERRM(9
0776      &           ,'(PYSLHA:) Missing A_b trilinear coupling')
0777             IF (MSPC(11).EQ.0) CALL PYERRM(9
0778      &           ,'(PYSLHA:) Missing A_t trilinear coupling')
0779             IF (MSPC(12).EQ.0) CALL PYERRM(9
0780      &           ,'(PYSLHA:) Missing A_tau trilinear coupling')
0781             IF (MSPC(17).LT.1) CALL PYERRM(9
0782      &           ,'(PYSLHA:) Missing Higgs mixing angle alpha')
0783           ENDIF
0784 C...Check wavefunction normalizations.
0785 C...Sfermions
0786           DO 320 ISPC=5,7
0787             IF (MSPC(ISPC).EQ.4) THEN
0788               KFSM=ISPC
0789               IF (ISPC.EQ.7) KFSM=15
0790               CHECK=ABS(SFMIX(KFSM,1)*SFMIX(KFSM,4)-SFMIX(KFSM,2)
0791      &             *SFMIX(KFSM,3))
0792               IF (ABS(1D0-CHECK).GT.1D-3) THEN
0793                 KCSM=PYCOMP(KFSM)
0794                 CALL PYERRM(17
0795      &               ,'(PYSLHA:) Non-orthonormal mixing matrix for ~'
0796      &               //CHAF(KCSM,1))
0797               ENDIF
0798             ENDIF
0799   320     CONTINUE
0800 C...Neutralinos + charginos
0801           DO 340 J=1,4
0802             CN1=0D0
0803             CN2=0D0
0804             CU1=0D0
0805             CU2=0D0
0806             CV1=0D0
0807             CV2=0D0
0808             DO 330 L=1,4
0809               CN1=CN1+ZMIX(J,L)**2
0810               CN2=CN2+ZMIX(L,J)**2
0811               IF (J.LE.2.AND.L.LE.2) THEN
0812                 CU1=CU1+UMIX(J,L)**2
0813                 CU2=CU2+UMIX(L,J)**2
0814                 CV1=CV1+VMIX(J,L)**2
0815                 CV2=CV2+VMIX(L,J)**2
0816               ENDIF
0817   330       CONTINUE
0818 C...NMIX normalization
0819             IF (MSPC(2).EQ.16.AND.(ABS(1D0-CN1).GT.1D-3.OR.ABS(1D0-CN2)
0820      &           .GT.1D-3).AND.IMSS(13).EQ.0) THEN
0821               CALL PYERRM(19,
0822      &             '(PYSLHA:) NMIX: Inconsistent normalization.')
0823               WRITE(MSTU(11),'(7x,I2,1x,":",2(1x,F7.4))') J, CN1, CN2
0824             ENDIF
0825 C...UMIX, VMIX normalizations
0826             IF (MSPC(3).EQ.4.OR.MSPC(4).EQ.4.AND.IMSS(13).EQ.0) THEN
0827               IF (J.LE.2) THEN
0828                 IF (ABS(1D0-CU1).GT.1D-3.OR.ABS(1D0-CU2).GT.1D-3) THEN
0829                   CALL PYERRM(19
0830      &                ,'(PYSLHA:) UMIX: Inconsistent normalization.')
0831                   WRITE(MSTU(11),'(7x,I2,1x,":",2(1x,F6.2))') J, CU1,
0832      &                 CU2
0833                 ENDIF
0834                 IF (ABS(1D0-CV1).GT.1D-3.OR.ABS(1D0-CV2).GT.1D-3) THEN
0835                   CALL PYERRM(19,
0836      &                '(PYSLHA:) VMIX: Inconsistent normalization.')
0837                   WRITE(MSTU(11),'(7x,I2,1x,":",2(1x,F6.2))') J, CV1,
0838      &                 CV2
0839                 ENDIF
0840               ENDIF
0841             ENDIF
0842   340     CONTINUE
0843           IF (MSTU(27).EQ.MSTU27.AND.MSTU(23).EQ.MSTU23) THEN
0844             WRITE(MSTU(11),'(1x,"*"/1x,A/1x,"*")')
0845      &           '*  PYSLHA:  No spectrum inconsistencies were found.'
0846           ELSE
0847             WRITE(MSTU(11),'(1x,"*"/1x,A/1x,"*",A/1x,"*",A/)')
0848      &           '* (PYSLHA:) INCONSISTENT SPECTRUM WARNING.'
0849      &           ,'Warning: one or more (serious)'//
0850      &           ' inconsistencies were found in the spectrum!!!'
0851      &           ,'Read the error messages above and check your'//
0852      &           ' input file.'
0853           ENDIF
0854 C...Increase precision in Higgs sector using FeynHiggs
0855           IF (IMSS(4).EQ.3) THEN
0856 C...FeynHiggs needs MSOFT.
0857             IERR=0
0858             IF (MSPC(18).EQ.0) THEN
0859               WRITE(MSTU(11),'(1x,"*"/1x,A/)')
0860      &             '* (PYSLHA:) BLOCK MSOFT not found in SLHA file.'//
0861      &              ' Cannot call FeynHiggs.'
0862               IERR=-1
0863             ELSE
0864               WRITE(MSTU(11),'(1x,/1x,A/)')
0865      &             '* (PYSLHA:) Now calling FeynHiggs.'
0866               CALL PYFEYN(IERR)
0867               IF (IERR.NE.0) IMSS(4)=2
0868             ENDIF
0869           ENDIF
0870         ELSEIF (MUPDA.EQ.2.AND.IRETRN.EQ.0) THEN
0871           KF=KFORIG
0872           KC=PYCOMP(KF)
0873           WRITE(CHKF,8300) KF
0874           IF(MIN(PMAS(KC,1),PMAS(KC,2),PMAS(KC,3),PMAS(KC,1)-PMAS(KC,3
0875      $         ),PMAS(KC,4)).LT.0D0.OR.MDCY(KC,3).LT.0.OR.(MDCY(KC,3)
0876      $         .EQ.0.AND.MDCY(KC,1).GE.1)) CALL PYERRM(17
0877      $         ,'(PYSLHA:) Mass/width/life/(# channels) wrong for KF='
0878      $         //CHKF)
0879           BRSUM=0D0
0880           BROPN=0D0
0881           DO 360 IDA=MDCY(KC,2),MDCY(KC,2)+MDCY(KC,3)-1
0882             IF(MDME(IDA,2).GT.80) GOTO 360
0883             KQ=KCHG(KC,1)
0884             PMS=PMAS(KC,1)-PMAS(KC,3)-PARJ(64)
0885             MERR=0
0886             DO 350 J=1,5
0887               KP=KFDP(IDA,J)
0888               IF(KP.EQ.0.OR.KP.EQ.81.OR.IABS(KP).EQ.82) THEN
0889                 IF(KP.EQ.81) KQ=0
0890               ELSEIF(PYCOMP(KP).EQ.0) THEN
0891                 MERR=3
0892               ELSE
0893                 KQ=KQ-PYCHGE(KP)
0894                 KPC=PYCOMP(KP)
0895                 PMS=PMS-PMAS(KPC,1)
0896                 IF(MSTJ(24).GT.0) PMS=PMS+0.5D0*MIN(PMAS(KPC,2),
0897      &               PMAS(KPC,3))
0898               ENDIF
0899   350       CONTINUE
0900             IF(KQ.NE.0) MERR=MAX(2,MERR)
0901             IF(MWID(KC).EQ.0.AND.KF.NE.311.AND.PMS.LT.0D0)
0902      &           MERR=MAX(1,MERR)
0903             IF(MERR.EQ.3) CALL PYERRM(17,
0904      &           '(PYSLHA:) Unknown particle code in decay of KF ='
0905      $           //CHKF)
0906             IF(MERR.EQ.2) CALL PYERRM(17,
0907      &           '(PYSLHA:) Charge not conserved in decay of KF ='
0908      $           //CHKF)
0909             IF(MERR.EQ.1) CALL PYERRM(7,
0910      &           '(PYSLHA:) Kinematically unallowed decay of KF ='
0911      $           //CHKF)
0912             BRSUM=BRSUM+BRAT(IDA)
0913             IF (MDME(IDA,1).GT.0) BROPN=BROPN+BRAT(IDA)
0914   360     CONTINUE
0915 C...Check branching ratio sum.
0916           IF (BROPN.LE.0D0) THEN
0917 C...If zero, set stable. 
0918              WRITE(CHTMP,8500) BROPN
0919              CALL PYERRM(7
0920      &            ,"(PYSLHA:) Effective BR sum for KF="//CHKF//' is '//
0921      &            CHTMP(9:16)//'. Changed to stable.')
0922              PMAS(KC,2)=1D-6
0923              MWID(KC)=0
0924 C...If BR's > 1, rescale.
0925           ELSEIF (BRSUM.GT.(1D0+1D-6)) THEN
0926              WRITE(CHTMP,8500) BRSUM
0927              CALL PYERRM(7
0928      &            ,"(PYSLHA:) Forced rescaling of BR's for KF="//CHKF//
0929      &            ' ; sum was'//CHTMP(9:16)//'.')
0930              FAC=1D0/BRSUM
0931              DO 370 IDA=MDCY(KC,2),MDCY(KC,2)+MDCY(KC,3)-1
0932                 IF(MDME(IDA,2).GT.80) GOTO 370
0933                 BRAT(IDA)=FAC*BRAT(IDA)
0934  370         CONTINUE
0935           ELSEIF (BRSUM.LT.(1D0-1D-6)) THEN
0936 C...If BR's < 1, insert dummy mode for proper cross section rescaling.
0937              WRITE(CHTMP,8500) BRSUM
0938              CALL PYERRM(7
0939      &            ,"(PYSLHA:) Sum of BR's for KF="//CHKF//' is '//
0940      &            CHTMP(9:16)//'. Dummy mode will be inserted.')
0941 C...  Insert dummy mode
0942              MDCY(KC,3)=MDCY(KC,3)+1
0943              IDA=MDCY(KC,2)+MDCY(KC,3)-1
0944              BRAT(IDA)=1D0-BRSUM
0945              KFDP(IDA,1)=0
0946              KFDP(IDA,2)=0
0947              KFDP(IDA,3)=0
0948              KFDP(IDA,4)=0
0949              KFDP(IDA,5)=0
0950              MDME(IDA,1)=0
0951              BRSUM=1D0
0952           ENDIF
0953        ENDIF
0954  
0955 C...WRITE SPECTRUM ON SLHA FILE
0956       ELSEIF(MUPDA.EQ.3) THEN
0957 C...If SPYTHIA or ISASUSY runtime was called for SUGRA, update PARMIN.
0958         IF (IMSS(1).EQ.2.OR.IMSS(1).EQ.12) THEN
0959           MODSEL(1)=1
0960           PARMIN(1)=RMSS(8)
0961           PARMIN(2)=RMSS(1)
0962           PARMIN(3)=RMSS(5)
0963           PARMIN(4)=SIGN(1D0,RMSS(4))
0964           PARMIN(5)=RMSS(36)
0965         ENDIF
0966 C...Write spectrum
0967         WRITE(LFN,7000) 'SLHA MSSM spectrum'
0968         WRITE(LFN,7000) 'Pythia 6.4: T. Sjostrand, S. Mrenna,'
0969      &    // ' P. Skands.'
0970         WRITE(LFN,7010) 'MODSEL',  'Model selection'
0971         WRITE(LFN,7110) 1, MODSEL(1)
0972         WRITE(LFN,7010) 'MINPAR', 'Parameters for minimal model.'
0973         IF (MODSEL(1).EQ.1) THEN
0974           WRITE(LFN,7210) 1, PARMIN(1), 'm0'
0975           WRITE(LFN,7210) 2, PARMIN(2), 'm12'
0976           WRITE(LFN,7210) 3, PARMIN(3), 'tan(beta)'
0977           WRITE(LFN,7210) 4, PARMIN(4), 'sign(mu)'
0978           WRITE(LFN,7210) 5, PARMIN(5), 'a0'
0979         ELSEIF(MODSEL(2).EQ.2) THEN
0980           WRITE(LFN,7210) 1, PARMIN(1), 'Lambda'
0981           WRITE(LFN,7210) 2, PARMIN(2), 'M'
0982           WRITE(LFN,7210) 3, PARMIN(3), 'tan(beta)'
0983           WRITE(LFN,7210) 4, PARMIN(4), 'sign(mu)'
0984           WRITE(LFN,7210) 5, PARMIN(5), 'N5'
0985           WRITE(LFN,7210) 6, PARMIN(6), 'c_grav'
0986         ENDIF
0987         WRITE(LFN,7000) ' '
0988         WRITE(LFN,7010) 'MASS', 'Mass spectrum'
0989         DO 380 I=1,36
0990           KF=KFSUSY(I)
0991           KC=PYCOMP(KF)
0992           IF (KF.EQ.1000039.AND.MODSEL(1).NE.2) GOTO 380
0993           KFSM=KF-KSUSY1
0994           IF (KFSM.GE.22.AND.KFSM.LE.37) THEN
0995             IF (KFSM.EQ.22)  WRITE(LFN,7220) KF, SMZ(1), CHAF(KC,1)
0996             IF (KFSM.EQ.23)  WRITE(LFN,7220) KF, SMZ(2), CHAF(KC,1)
0997             IF (KFSM.EQ.25)  WRITE(LFN,7220) KF, SMZ(3), CHAF(KC,1)
0998             IF (KFSM.EQ.35)  WRITE(LFN,7220) KF, SMZ(4), CHAF(KC,1)
0999             IF (KFSM.EQ.24)  WRITE(LFN,7220) KF, SMW(1), CHAF(KC,1)
1000             IF (KFSM.EQ.37)  WRITE(LFN,7220) KF, SMW(2), CHAF(KC,1)
1001           ELSE
1002             WRITE(LFN,7220) KF, PMAS(KC,1), CHAF(KC,1)
1003           ENDIF
1004   380   CONTINUE
1005 C...SUSY scale
1006         RMSUSY=SQRT(PMAS(PYCOMP(KSUSY1+6),1)*PMAS(PYCOMP(KSUSY2+6),1))
1007         WRITE(LFN,7020) 'HMIX',RMSUSY,'Higgs parameters'
1008         WRITE(LFN,7210) 1, RMSS(4),'mu'
1009         WRITE(LFN,7010) 'ALPHA',' '
1010         WRITE(LFN,7210) 1, RMSS(18), 'alpha'
1011         WRITE(LFN,7020) 'AU',RMSUSY
1012         WRITE(LFN,7410) 3, 3, RMSS(16), 'A_t'
1013         WRITE(LFN,7020) 'AD',RMSUSY
1014         WRITE(LFN,7410) 3, 3, RMSS(15), 'A_b'
1015         WRITE(LFN,7020) 'AE',RMSUSY
1016         WRITE(LFN,7410) 3, 3, RMSS(17), 'A_tau'
1017         WRITE(LFN,7010) 'STOPMIX','~t mixing matrix'
1018         WRITE(LFN,7410) 1, 1, SFMIX(6,1)
1019         WRITE(LFN,7410) 1, 2, SFMIX(6,2)
1020         WRITE(LFN,7410) 2, 1, SFMIX(6,3)
1021         WRITE(LFN,7410) 2, 2, SFMIX(6,4)
1022         WRITE(LFN,7010) 'SBOTMIX','~b mixing matrix'
1023         WRITE(LFN,7410) 1, 1, SFMIX(5,1)
1024         WRITE(LFN,7410) 1, 2, SFMIX(5,2)
1025         WRITE(LFN,7410) 2, 1, SFMIX(5,3)
1026         WRITE(LFN,7410) 2, 2, SFMIX(5,4)
1027         WRITE(LFN,7010) 'STAUMIX','~tau mixing matrix'
1028         WRITE(LFN,7410) 1, 1, SFMIX(15,1)
1029         WRITE(LFN,7410) 1, 2, SFMIX(15,2)
1030         WRITE(LFN,7410) 2, 1, SFMIX(15,3)
1031         WRITE(LFN,7410) 2, 2, SFMIX(15,4)
1032         WRITE(LFN,7010) 'NMIX','~chi0 mixing matrix'
1033         DO 400 I1=1,4
1034           DO 390 I2=1,4
1035             WRITE(LFN,7410) I1, I2, ZMIX(I1,I2)
1036   390     CONTINUE
1037   400   CONTINUE
1038         WRITE(LFN,7010) 'UMIX','~chi^+ U mixing matrix'
1039         DO 420 I1=1,2
1040           DO 410 I2=1,2
1041             WRITE(LFN,7410) I1, I2, UMIX(I1,I2)
1042   410     CONTINUE
1043   420   CONTINUE
1044         WRITE(LFN,7010) 'VMIX','~chi^+ V mixing matrix'
1045         DO 440 I1=1,2
1046           DO 430 I2=1,2
1047             WRITE(LFN,7410) I1, I2, VMIX(I1,I2)
1048   430     CONTINUE
1049   440   CONTINUE
1050         WRITE(LFN,7010) 'SPINFO'
1051         IF (IMSS(1).EQ.2) THEN
1052           CPRO(1)='PYTHIA'
1053           CVER(1)='6.4'
1054         ELSEIF (IMSS(1).EQ.12) THEN
1055           ISAVER=VISAJE()
1056           CPRO(1)='ISASUSY'
1057           CVER(1)=ISAVER(1:12)
1058         ENDIF
1059         WRITE(LFN,7310) 1, CPRO(1), 'Spectrum Calculator'
1060         WRITE(LFN,7310) 2, CVER(1), 'Version number'
1061       ENDIF
1062  
1063 C...Print user information about spectrum
1064       IF (MUPDA.EQ.1.OR.MUPDA.EQ.3) THEN
1065         IF (CPRO(MOD(MUPDA,2)).NE.' '.AND.CVER(MOD(MUPDA,2)).NE.' ')
1066      &       WRITE(MSTU(11),5030) CPRO(1), CVER(1)
1067         IF (IMSS(4).EQ.3) WRITE(MSTU(11),5040)
1068         IF (MUPDA.EQ.1) THEN
1069           WRITE(MSTU(11),5020) LFN
1070         ELSE
1071           WRITE(MSTU(11),5010) LFN
1072         ENDIF
1073  
1074         WRITE(MSTU(11),5400)
1075         WRITE(MSTU(11),5500) 'Pole masses'
1076         WRITE(MSTU(11),5700) (RMFUN(KSUSY1+IP),IP=1,6)
1077      $       ,(RMFUN(KSUSY2+IP),IP=1,6)
1078         WRITE(MSTU(11),5800) (RMFUN(KSUSY1+IP),IP=11,16)
1079      $       ,(RMFUN(KSUSY2+IP),IP=11,16)
1080         IF (IMSS(13).EQ.0) THEN
1081           WRITE(MSTU(11),5900) RMFUN(KSUSY1+21),RMFUN(KSUSY1+22)
1082      $         ,RMFUN(KSUSY1+23),RMFUN(KSUSY1+25),RMFUN(KSUSY1+35),
1083      $         RMFUN(KSUSY1+24),RMFUN(KSUSY1+37)
1084           WRITE(MSTU(11),6000) CHAF(25,1),CHAF(35,1),CHAF(36,1),
1085      &         CHAF(37,1), ' ', ' ',' ',' ',
1086      &         RMFUN(25), RMFUN(35), RMFUN(36), RMFUN(37)
1087         ELSEIF (IMSS(13).EQ.1) THEN
1088           KF1=KSUSY1+21
1089           KF2=KSUSY1+22
1090           KF3=KSUSY1+23
1091           KF4=KSUSY1+25
1092           KF5=KSUSY1+35
1093           KF6=KSUSY1+45
1094           KF7=KSUSY1+24
1095           KF8=KSUSY1+37
1096           WRITE(MSTU(11),6000) CHAF(PYCOMP(KF1),1),CHAF(PYCOMP(KF2),1),
1097      &         CHAF(PYCOMP(KF3),1),CHAF(PYCOMP(KF4),1),
1098      &         CHAF(PYCOMP(KF5),1),CHAF(PYCOMP(KF6),1),
1099      &         CHAF(PYCOMP(KF7),1),CHAF(PYCOMP(KF8),1),
1100      &         RMFUN(KF1),RMFUN(KF2),RMFUN(KF3),RMFUN(KF4),
1101      &         RMFUN(KF5),RMFUN(KF6),RMFUN(KF7),RMFUN(KF8)
1102           WRITE(MSTU(11),6000) CHAF(25,1), CHAF(35,1), CHAF(45,1),
1103      &         CHAF(36,1), CHAF(46,1), CHAF(37,1),' ',' ',
1104      &         RMFUN(25), RMFUN(35), RMFUN(45), RMFUN(36), RMFUN(46),
1105      &         RMFUN(37)
1106         ENDIF
1107         WRITE(MSTU(11),5400)
1108         WRITE(MSTU(11),5500) 'Mixing structure'
1109         WRITE(MSTU(11),6100) ((ZMIX(I,J), J=1,4),I=1,4)
1110         WRITE(MSTU(11),6200) (UMIX(1,J), J=1,2),(VMIX(1,J),J=1,2)
1111      &       ,(UMIX(2,J), J=1,2),(VMIX(2,J),J=1,2)
1112         WRITE(MSTU(11),6300) (SFMIX(5,J), J=1,2),(SFMIX(6,J),J=1,2)
1113      &       ,(SFMIX(15,J), J=1,2),(SFMIX(5,J),J=3,4),(SFMIX(6,J), J=3,4
1114      &       ),(SFMIX(15,J),J=3,4)
1115         WRITE(MSTU(11),5400)
1116         WRITE(MSTU(11),5500) 'Couplings'
1117         WRITE(MSTU(11),6400) RMSS(15),RMSS(16),RMSS(17)
1118         WRITE(MSTU(11),6450) RMSS(18), RMSS(5), RMSS(4)
1119         WRITE(MSTU(11),5400)
1120         WRITE(MSTU(11),6500)
1121  
1122       ENDIF
1123  
1124 C...Only rewind when reading
1125       IF (MUPDA.LE.2.OR.MUPDA.EQ.5) REWIND(LFN)
1126  
1127  9999 RETURN
1128  
1129 C...Serious error catching
1130   460 write(*,*) '* (PYSLHA:) read BLOCK error on line',NLINE
1131       write(*,*) CHINL(1:80)
1132       CALL PYSTOP(106)
1133   470 WRITE(*,*) '* (PYSLHA:) read DECAY error on line',NLINE
1134       WRITE(*,*) CHINL(1:72)
1135       CALL PYSTOP(106)
1136   480 WRITE(*,*) '* (PYSLHA:) read BR error on line',NLINE
1137       WRITE(*,*) CHINL(1:80)
1138       CALL PYSTOP(106)
1139   490 WRITE(*,*) '* (PYSLHA:) read NDA error on line',NLINE
1140       WRITE(*,*) CHINL(1:80)
1141       CALL PYSTOP(106)
1142   500 WRITE(*,*) '* (PYSLHA:) decay daughter read error on line',NLINE
1143       WRITE(*,*) CHINL(1:80)
1144   510 WRITE(*,*) '* (PYSLHA:) read Q error in BLOCK ',CHBLCK
1145       CALL PYSTOP(106)
1146   520 WRITE(*,*) '* (PYSLHA:) read error in line ',NLINE,':'
1147       WRITE(*,*) CHINL(1:80)
1148       CALL PYSTOP(106)
1149  
1150  8300 FORMAT(I9)
1151  8500 FORMAT(F16.5)
1152  
1153 C...Formats for user information printout.
1154  5000 FORMAT(1x,15('*'),1x,'PYSLHA v1.09: SUSY/BSM SPECTRUM '
1155      &     ,'INTERFACE',1x,15('*')/1x,'*',2x
1156      &     ,'PYSLHA:  Last Change',1x,A,1x,'-',1x,'P.Z. Skands')
1157  5010 FORMAT(1x,'*',3x,'Wrote spectrum file on unit: ',I3)
1158  5020 FORMAT(1x,'*',3x,'Read spectrum file on unit: ',I3)
1159  5030 FORMAT(1x,'*',3x,'Spectrum Calculator was: ',A,' version ',A)
1160  5040 FORMAT(1x,'*',3x,'Higgs sector corrected with FeynHiggs')
1161  5100 FORMAT(1x,'*',1x,'Model parameters:'/1x,'*',1x,'----------------')
1162  5200 FORMAT(1x,'*',1x,3x,'M_0',6x,'M_1/2',5x,'A_0',3x,'Tan(beta)',
1163      &     3x,'Sgn(mu)',3x,'M_t'/1x,'*',1x,4(F8.2,1x),I8,2x,F8.2)
1164  5300 FORMAT(1x,'*'/1x,'*',1x,'Model spectrum :'/1x,'*',1x
1165      &     ,'----------------')
1166  5400 FORMAT(1x,'*',1x,A)
1167  5500 FORMAT(1x,'*',1x,A,':')
1168  5600 FORMAT(1x,'*',2x,2x,'M_GUT',2x,2x,'g_GUT',2x,1x,'alpha_GUT'/
1169      &       1x,'*',2x,1P,2(1x,E8.2),2x,E8.2)
1170  5700 FORMAT(1x,'*',4x,4x,'~d',2x,1x,4x,'~u',2x,1x,4x,'~s',2x,1x,
1171      &     4x,'~c',2x,1x,1x,'~b(12)',1x,1x,1x,'~t(12)'/1x,'*',2x,'L',1x
1172      &     ,6(F8.2,1x)/1x,'*',2x,'R',1x,6(F8.2,1x))
1173  5800 FORMAT(1x,'*'/1x,'*',4x,4x,'~e',2x,1x,2x,'~nu_e',2x,1x,3x,'~mu',2x
1174      &     ,1x,1x,'~nu_mu',1x,1x,'~tau(12)',1x,1x,'~nu_tau'/1x,'*',2x
1175      &     ,'L',1x,6(F8.2,1x)/1x,'*',2x,'R',1x,6(F8.2,1x))
1176  5900 FORMAT(1x,'*'/1x,'*',4x,4x,'~g',2x,1x,1x,'~chi_10',1x,1x,'~chi_20'
1177      &     ,1x,1x,'~chi_30',1x,1x,'~chi_40',1x,1x,'~chi_1+',1x
1178      &     ,1x,'~chi_2+'/1x,'*',3x,1x,7(F8.2,1x))
1179  6000 FORMAT(1x,'*'/1x,'*',3x,1x,8(1x,A7,1x)/1x,'*',3x,1x,8(F8.2,1x))
1180  6100 FORMAT(1x,'*',11x,'|',3x,'~B',3x,'|',2x,'~W_3',2x,'|',2x
1181      &     ,'~H_1',2x,'|',2x,'~H_2',2x,'|'/1x,'*',3x,'~chi_10',1x,4('|'
1182      &     ,1x,F6.3,1x),'|'/1x,'*',3x,'~chi_20',1x,4('|'
1183      &     ,1x,F6.3,1x),'|'/1x,'*',3x,'~chi_30',1x,4('|'
1184      &     ,1x,F6.3,1x),'|'/1x,'*',3x,'~chi_40',1x,4('|'
1185      &     ,1x,F6.3,1x),'|')
1186  6200 FORMAT(1x,'*'/1x,'*',6x,'L',4x,'|',3x,'~W',3x,'|',3x,'~H',3x,'|'
1187      &     ,12x,'R',4x,'|',3x,'~W',3x,'|',3x,'~H',3x,'|'/1x,'*',3x
1188      &     ,'~chi_1+',1x,2('|',1x,F6.3,1x),'|',9x,'~chi_1+',1x,2('|',1x
1189      &     ,F6.3,1x),'|'/1x,'*',3x,'~chi_2+',1x,2('|',1x,F6.3,1x),'|',9x
1190      &     ,'~chi_2+',1x,2('|',1x,F6.3,1x),'|')
1191  6300 FORMAT(1x,'*'/1x,'*',8x,'|',2x,'~b_L',2x,'|',2x,'~b_R',2x,'|',8x
1192      &     ,'|',2x,'~t_L',2x,'|',2x,'~t_R',2x,'|',10x
1193      &     ,'|',1x,'~tau_L',1x,'|',1x,'~tau_R',1x,'|'/
1194      &     1x,'*',3x,'~b_1',1x,2('|',1x,F6.3,1x),'|',3x,'~t_1',1x,2('|'
1195      &     ,1x,F6.3,1x),'|',3x,'~tau_1',1x,2('|',1x,F6.3,1x),'|'/
1196      &     1x,'*',3x,'~b_2',1x,2('|',1x,F6.3,1x),'|',3x,'~t_2',1x,2('|'
1197      &     ,1x,F6.3,1x),'|',3x,'~tau_2',1x,2('|',1x,F6.3,1x),'|')
1198  6400 FORMAT(1x,'*',3x,'  A_b = ',F8.2,4x,'      A_t = ',F8.2,4x
1199      &     ,'A_tau = ',F8.2)
1200  6450 FORMAT(1x,'*',3x,'alpha = ',F8.2,4x,'tan(beta) = ',F8.2,4x
1201      &     ,'   mu = ',F8.2)
1202  6500 FORMAT(1x,32('*'),1x,'END OF PYSLHA',1x,31('*'))
1203  
1204 C...Format to use for comments
1205  7000 FORMAT('# ',A)
1206 C...Format to use for block statements
1207  7010 FORMAT('Block',1x,A,3x,'#',1x,A)
1208  7020 FORMAT('Block',1x,A,1x,'Q=',1P,E16.8,0P,3x,'#',1x,A)
1209 C...Indexed Int
1210  7110 FORMAT(1x,I4,1x,I4,3x,'#')
1211 C...Non-Indexed Double
1212  7200 FORMAT(9x,1P,E16.8,0P,3x,'#',1x,A)
1213 C...Indexed Double
1214  7210 FORMAT(1x,I4,3x,1P,E16.8,0P,3x,'#',1x,A)
1215 C...Long Indexed Double (PDG + double)
1216  7220 FORMAT(1x,I9,3x,1P,E16.8,0P,3x,'#',1x,A)
1217 C...Indexed Char(12)
1218  7310 FORMAT(1x,I4,3x,A12,3x,'#',1x,A)
1219 C...Single matrix
1220  7410 FORMAT(1x,I2,1x,I2,3x,1P,E16.8,0P,3x,'#',1x,A)
1221 C...Double Matrix
1222  7420 FORMAT(1x,I2,1x,I2,3x,1P,E16.8,3x,E16.8,0P,3x,'#',1x,A)
1223 C...Write Decay Table
1224  7500 FORMAT('Decay',1x,I9,1x,'WIDTH=',1P,E16.8,0P,3x,'#',1x,A)
1225  7510 FORMAT(4x,1P,E16.8,0P,3x,I2,3x,'IDA=',1x,5(1x,I9),3x,'#',1x,A)
1226  
1227       END
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!