| [2403] | 1 | * | 
|---|
|  | 2 | * $Id: minuit.F,v 1.1.1.1 2003-06-11 14:18:26 cmv Exp $ | 
|---|
|  | 3 | * | 
|---|
|  | 4 | * $Log: not supported by cvs2svn $ | 
|---|
|  | 5 | * Revision 1.1.1.1  1996/03/07 14:31:28  mclareni | 
|---|
|  | 6 | * Minuit | 
|---|
|  | 7 | * | 
|---|
|  | 8 | * | 
|---|
|  | 9 | #include "minuit/pilot.h" | 
|---|
|  | 10 | SUBROUTINE MINUIT(FCN,FUTIL) | 
|---|
|  | 11 | #include "minuit/d506dp.inc" | 
|---|
|  | 12 | #include "minuit/d506cm.inc" | 
|---|
|  | 13 | C | 
|---|
|  | 14 | C  CPNAM   Parameter name (10 characters) | 
|---|
|  | 15 | C  U       External (visible to user in FCN) value of parameter | 
|---|
|  | 16 | C  ALIM, BLIM Lower and upper parameter limits. If both zero, no limits. | 
|---|
|  | 17 | C  ERP,ERN Positive and negative MINOS errors, if calculated. | 
|---|
|  | 18 | C  WERR    External parameter error (standard deviation, defined by UP) | 
|---|
|  | 19 | C  GLOBCC  Global Correlation Coefficient | 
|---|
|  | 20 | C  NVARL   =-1 if parameter undefined,      =0 if constant, | 
|---|
|  | 21 | C          = 1 if variable without limits,  =4 if variable with limits | 
|---|
|  | 22 | C   (Note that if parameter has been fixed, NVARL=1 or =4, and NIOFEX=0) | 
|---|
|  | 23 | C  NIOFEX  Internal parameter number, or zero if not currently variable | 
|---|
|  | 24 | C  NEXOFI  External parameter number for currently variable parameters | 
|---|
|  | 25 | C  X, XT   Internal parameter values (X are sometimes saved in XT) | 
|---|
|  | 26 | C  DIRIN   (Internal) step sizes for current step | 
|---|
|  | 27 | C  variables with names ending in ..S are saved values for fixed params | 
|---|
|  | 28 | C  VHMAT   (Internal) error matrix stored as Half MATrix, since | 
|---|
|  | 29 | C                it is symmetric | 
|---|
|  | 30 | C  VTHMAT  VHMAT is sometimes saved in VTHMAT, especially in MNMNOT | 
|---|
|  | 31 | C | 
|---|
|  | 32 | C  ISW definitions: | 
|---|
|  | 33 | C      ISW(1) =0 normally, =1 means CALL LIMIT EXCEEDED | 
|---|
|  | 34 | C      ISW(2) =0 means no error matrix | 
|---|
|  | 35 | C             =1 means only approximate error matrix | 
|---|
|  | 36 | C             =2 means full error matrix, but forced pos-def. | 
|---|
|  | 37 | C             =3 means good normal full error matrix exists | 
|---|
|  | 38 | C      ISW(3) =0 if Minuit is calculating the first derivatives | 
|---|
|  | 39 | C             =1 if first derivatives calculated inside FCN | 
|---|
|  | 40 | C      ISW(4) =-1 if most recent minimization did not converge. | 
|---|
|  | 41 | C             = 0 if problem redefined since most recent minimization. | 
|---|
|  | 42 | C             =+1 if most recent minimization did converge. | 
|---|
|  | 43 | C      ISW(5) is the PRInt level.  See SHO PRIntlevel | 
|---|
|  | 44 | C      ISW(6) = 0 for batch mode, =1 for interactive mode | 
|---|
|  | 45 | C                      =-1 for originally interactive temporarily batch | 
|---|
|  | 46 | C | 
|---|
|  | 47 | C  LWARN is true if warning messges are to be put out (default=true) | 
|---|
|  | 48 | C            SET WARN turns it on, set NOWarn turns it off | 
|---|
|  | 49 | C  LREPOR is true if exceptional conditions are put out (default=false) | 
|---|
|  | 50 | C            SET DEBUG turns it on, SET NODebug turns it off | 
|---|
|  | 51 | C  LIMSET is true if a parameter is up against limits (for MINOS) | 
|---|
|  | 52 | C  LNOLIM is true if there are no limits on any parameters (not yet used) | 
|---|
|  | 53 | C  LNEWMN is true if the previous process has unexpectedly improved FCN | 
|---|
|  | 54 | C  LPHEAD is true if a heading should be put out for the next parameter | 
|---|
|  | 55 | C        definition, false if a parameter has just been defined | 
|---|
|  | 56 | C | 
|---|
|  | 57 | EXTERNAL FCN,FUTIL | 
|---|
|  | 58 | CHARACTER*40 CWHYXT | 
|---|
|  | 59 | DATA CWHYXT/'FOR UNKNOWN REASONS                     '/ | 
|---|
|  | 60 | DATA JSYSRD,JSYSWR,JSYSSA/5,6,7/ | 
|---|
|  | 61 | C                                 . . . . . . . . . . initialize minuit | 
|---|
|  | 62 | WRITE (JSYSWR,'(1X,75(1H*))') | 
|---|
|  | 63 | CALL MNINIT (JSYSRD,JSYSWR,JSYSSA) | 
|---|
|  | 64 | C                                      . . . . initialize new data block | 
|---|
|  | 65 | 100 CONTINUE | 
|---|
|  | 66 | WRITE (ISYSWR,'(1X,75(1H*))') | 
|---|
|  | 67 | NBLOCK = NBLOCK + 1 | 
|---|
|  | 68 | WRITE (ISYSWR,'(26X,A,I4)')  'MINUIT DATA BLOCK NO.',NBLOCK | 
|---|
|  | 69 | WRITE (ISYSWR,'(1X,75(1H*))') | 
|---|
|  | 70 | C               . . . . . . . . . . .   set parameter lists to undefined | 
|---|
|  | 71 | CALL MNCLER | 
|---|
|  | 72 | C                                             . . . . . . . . read title | 
|---|
|  | 73 | CALL MNREAD(FCN,1,IFLGUT,FUTIL) | 
|---|
|  | 74 | IF (IFLGUT .EQ. 2)  GO TO 500 | 
|---|
|  | 75 | IF (IFLGUT .EQ. 3)  GO TO 600 | 
|---|
|  | 76 | C                                        . . . . . . . . read parameters | 
|---|
|  | 77 | CALL MNREAD(FCN,2,IFLGUT,FUTIL) | 
|---|
|  | 78 | IF (IFLGUT .EQ. 2)  GO TO 500 | 
|---|
|  | 79 | IF (IFLGUT .EQ. 3)  GO TO 600 | 
|---|
|  | 80 | IF (IFLGUT .EQ. 4)  GO TO 700 | 
|---|
|  | 81 | C                              . . . . . . verify FCN not time-dependent | 
|---|
|  | 82 | WRITE (ISYSWR,'(/A,A)') ' MINUIT: FIRST CALL TO USER FUNCTION,', | 
|---|
|  | 83 | +    ' WITH IFLAG=1' | 
|---|
|  | 84 | NPARX = NPAR | 
|---|
|  | 85 | CALL MNINEX(X) | 
|---|
|  | 86 | FZERO = UNDEFI | 
|---|
|  | 87 | CALL FCN(NPARX,GIN,FZERO,U,1,FUTIL) | 
|---|
|  | 88 | FIRST = UNDEFI | 
|---|
|  | 89 | CALL FCN(NPARX,GIN,FIRST,U,4,FUTIL) | 
|---|
|  | 90 | NFCN = 2 | 
|---|
|  | 91 | IF (FZERO.EQ.UNDEFI .AND. FIRST.EQ.UNDEFI)  THEN | 
|---|
|  | 92 | CWHYXT = 'BY ERROR IN USER FUNCTION.  ' | 
|---|
|  | 93 | WRITE (ISYSWR,'(/A,A/)') ' USER HAS NOT CALCULATED FUNCTION', | 
|---|
|  | 94 | +    ' VALUE WHEN IFLAG=1 OR 4' | 
|---|
|  | 95 | GO TO 800 | 
|---|
|  | 96 | ENDIF | 
|---|
|  | 97 | AMIN = FIRST | 
|---|
|  | 98 | IF (FIRST .EQ. UNDEFI) AMIN=FZERO | 
|---|
|  | 99 | CALL MNPRIN(1,AMIN) | 
|---|
|  | 100 | NFCN = 2 | 
|---|
|  | 101 | IF (FIRST .EQ. FZERO)  GO TO 300 | 
|---|
|  | 102 | FNEW = 0.0 | 
|---|
|  | 103 | CALL FCN(NPARX,GIN,FNEW,U,4,FUTIL) | 
|---|
|  | 104 | IF  (FNEW .NE. AMIN) WRITE (ISYSWR,280) AMIN, FNEW | 
|---|
|  | 105 | 280 FORMAT (/' MINUIT WARNING: PROBABLE ERROR IN USER FUNCTION.'/ | 
|---|
|  | 106 | +         ' FOR FIXED VALUES OF PARAMETERS, FCN IS TIME-DEPENDENT'/ | 
|---|
|  | 107 | +         ' F =',E22.14,' FOR FIRST CALL'/ | 
|---|
|  | 108 | +         ' F =',E22.14,' FOR SECOND CALL.'/) | 
|---|
|  | 109 | NFCN = 3 | 
|---|
|  | 110 | 300 FVAL3 = 2.0*AMIN+1.0 | 
|---|
|  | 111 | C                                   . . . . . . . . . . . read commands | 
|---|
|  | 112 | CALL MNREAD(FCN,3,IFLGUT,FUTIL) | 
|---|
|  | 113 | IF (IFLGUT .EQ. 2)  GO TO 500 | 
|---|
|  | 114 | IF (IFLGUT .EQ. 3)  GO TO 600 | 
|---|
|  | 115 | IF (IFLGUT .EQ. 4)  GO TO 700 | 
|---|
|  | 116 | CWHYXT = 'BY MINUIT COMMAND: '//CWORD | 
|---|
|  | 117 | IF (INDEX(CWORD,'STOP').GT. 0)  GO TO 800 | 
|---|
|  | 118 | IF (INDEX(CWORD,'EXI') .GT. 0)  GO TO 800 | 
|---|
|  | 119 | IF (INDEX(CWORD,'RET') .EQ. 0)  GO TO 100 | 
|---|
|  | 120 | CWHYXT = 'AND RETURNS TO USER PROGRAM.    ' | 
|---|
|  | 121 | WRITE (ISYSWR,'(A,A)')  ' ..........MINUIT TERMINATED ',CWHYXT | 
|---|
|  | 122 | RETURN | 
|---|
|  | 123 | C                                           . . . . . . stop conditions | 
|---|
|  | 124 | 500 CONTINUE | 
|---|
|  | 125 | CWHYXT = 'BY END-OF-DATA ON PRIMARY INPUT FILE.   ' | 
|---|
|  | 126 | GO TO 800 | 
|---|
|  | 127 | 600 CONTINUE | 
|---|
|  | 128 | CWHYXT = 'BY UNRECOVERABLE READ ERROR ON INPUT.   ' | 
|---|
|  | 129 | GO TO 800 | 
|---|
|  | 130 | 700 CONTINUE | 
|---|
|  | 131 | CWHYXT = ': FATAL ERROR IN PARAMETER DEFINITIONS. ' | 
|---|
|  | 132 | 800 WRITE (ISYSWR,'(A,A)')  ' ..........MINUIT TERMINATED ',CWHYXT | 
|---|
|  | 133 | STOP | 
|---|
|  | 134 | C | 
|---|
|  | 135 | C  ......................entry to set unit numbers  - - - - - - - - - - | 
|---|
|  | 136 | ENTRY MINTIO(I1,I2,I3) | 
|---|
|  | 137 | JSYSRD = I1 | 
|---|
|  | 138 | JSYSWR = I2 | 
|---|
|  | 139 | JSYSSA = I3 | 
|---|
|  | 140 | RETURN | 
|---|
|  | 141 | END | 
|---|