GLK.f

00001 *//////////////////////////////////////////////////////////////////////////////
00002 *//                                                                          //
00003 *//                       Pseudo-Class  GLK                                  //
00004 *//                                                                          //
00005 *//////////////////////////////////////////////////////////////////////////////
00006 *
00007 *
00008 *//////////////////////////////////////////////////////////////////////////////
00009 *// =======================================================================  //
00010 *// ==========================  _GLK_  ====================================  //
00011 *// ========== General Library of histogramming/ploting utilities =========  //
00012 *// ========== It is similar but not identical to HBOOK and HPLOT =========  //
00013 *// =======================================================================  //
00014 *//////////////////////////////////////////////////////////////////////////////
00015 *//                                                                          //
00016 *//                      Version:    1.30                                    //
00017 *//              Last correction:    January 1999                            //
00018 *//                                                                          //
00019 *//////////////////////////////////////////////////////////////////////////////
00020 *//
00021 *//    ********************************************************************
00022 *//    *  History of the package:                                         *
00023 *//    *  MINI-HBOOK writen by S. Jadach, Rutherford Lab. 1976            *
00024 *//    *  Rewritten December 1989 (S.J.)   and in 1997 (S.J.)             *
00025 *//    ********************************************************************
00026 *//
00027 *//  Installation remarks:
00028 *//  (1) printing backslash character depends on F77 compilator,
00029 *//      user may need to modify definition of BS variable in HPLCAP
00030 *//
00031 *//  Usage of the program:
00032 *//  (1) In many cases names and meanings of programs and their
00033 *//      parameters is the same as in original CERN libraries HBOOK
00034 *//  (2) Unlike to original HBOOK and HPLOT, all floating parameters
00035 *//      of the programs are in DOUBLE PRECISION !
00036 *//  (3) GLK stores histograms in DOUBLE PRECISION  and always with
00037 *//      errors. DOUBLE PRECISION  storage is essential for 10**7 events statistics!
00038 *//  (4) Output from GLK is a picture recorded as regular a LaTeX file
00039 *//      with frame and curves/histograms, it is easy to change fonts
00040 *//      add captions, merge plots, etc. by normal editing. Finally,
00041 *//      picture may be inserted in any place into LaTeX source of the
00042 *//      article.
00043 *//  (5) WARNING: two-dimensional histograms are not active!!!
00044 *//  
00045 *//////////////////////////////////////////////////////////////////////////////
00046 *//     List of procedures,  non-user subprograms in brackets                //
00047 *//////////////////////////////////////////////////////////////////////////////
00048 *    SUBR/FUNC       1 PAR. 2 PAR. 3 PAR. 4 PAR. 5 PAR. 6 PAR.
00049 *  ====================================================================
00050 *     (GLK_Initialize) ----   ----    ----   ----   ----   ----
00051 *      GLK_hi          INT    INT     ----   ----   ----   ----
00052 *      GLK_hie         INT    INT     ----   ----   ----   ----
00053 *      GLK_Fil1        INT    DBL     DBL    ----   ----   ----
00054 *      GLK_Fil2        INT    DBL     DBL    DBL    ----   ----
00055 *      GLK_Book1       INT    CHR*80  INT    DBL    DBL    ----
00056 *     (GLK_OptOut)     INT    INT     INT    INT    INT     INT
00057 * (L.F. GLK_Exist)     INT    -----  ------  ----   ----   ----
00058 *      GLK_Idopt       INT    CHR*4   -----  ----   ----   ----
00059 *      GLK_BookFun1    INT    CHR*80   INT   DBL    DBL  DP-FUNC
00060 *      GLK_Idopt       INT    CHR*4   -----  ----   ----   ----
00061 *      GLK_Book2       INT    CHR*80   INT   DBL    DBL     INT  DBL DBL
00062 *      GLK_PrintAll    ---    ----   ----   ----   ----   ----
00063 *      GLK_SetNout     INT    ----   ----   ----   ----   ----
00064 *      GLK_Print       INT    ----   ----   ----   ----   ----
00065 *      GLK_Operat      INT    CHR*1   INT    INT    DBL    DBL
00066 *      GLK_Hinbo1      INT    CHR*8   INT    DBL    DBL    ----
00067 *      GLK_Unpak       INT    DBL(*) CHR*(*) INT    ---    ----
00068 *      GLK_Pak         INT    DBL(*)  ----   ----   ---    ----
00069 *      GLK_Pake        INT    DBL(*)  ----   ----   ---    ----
00070 *      GLK_Range1      INT    DBL     DBL    ----   ---    ----
00071 *      GLK_Hinbo2      INT    INT     DBL    DBL    INT    DBL   DBL
00072 *      GLK_Ymaxim      INT    DBL     ----   ----   ---    ----
00073 *      GLK_Yminim      INT    DBL     ----   ----   ---    ----
00074 *      GLK_Reset       INT   CHR*(*)  ----   ----   ---    ----
00075 *      GLK_Delet       INT     ----   ----   ----   ----   ----
00076 *     (GLK_Copch)     CHR*80 CHR*80  ----   ----   ----   ----
00077 *     (GLK_hadres)     INT    INT   ----   ----   ----   ----
00078 *      GLK_Hrfile      INT   CHR*(*) CHR*(*) ----   ----   ----
00079 *      GLK_Hrout       INT    INT    CHR*8   ----   ----   ----
00080 *      GLK_Hrin        INT    INT     INT    ----   ----   ----
00081 *      GLK_Hrend     CHR*(*) ----    ----   ----   ----   ----
00082 *  *******************  HPLOT entries ******************
00083 *      GLK_PlInt       INT    ----    ----   ----   ----   ----
00084 *      GLK_PlCap       INT    ----    ----   ----   ----   ----
00085 *      GLK_PlEnd       ----   ----    ----   ----   ----   ----
00086 *      GLK_Plot        INT    CHR*1   CHR*1   INT   ----   ----
00087 *     (GLK_Plfram1)    INT      INT     INT  ----   ----   ----
00088 *     (GLK_SAxisX)     INT      DBL     DBL   INT    DBL   ----
00089 *     (GLK_SAxisY)     INT      DBL     DBL   INT    DBL   ----
00090 *     (GLK_PlHist)     INT      INT     DBL   DBL    INT    INT
00091 *     (GLK_PlHis2)     INT      INT     DBL   DBL    INT    INT
00092 *     (GLK_PlCirc)     INT      INT     INT   DBL    DBL    DBL
00093 *     (GLK_aprof)      DBL      INT     DBL  ----   ----   ----
00094 *      GLK_PlSet       INT      DBL    ----  ----   ----   ----
00095 *      GLK_PlTitle     INT    CHR*80   ----  ----   ----   ----
00096 *  *******************  WMONIT entries ******************
00097 *      GLK_WtMon       INT ???
00098 *  *******************************************************************
00099 *                         END OF TABLE
00100 *  *******************************************************************
00101 *          Map of memory for single histogram
00102 *          ----------------------------------
00103 *  (1-7) Header
00104 *  ist +1   mark      9999999999999
00105 *  ist +2   mark      9d12 + id*10 + 9
00106 *  ist +3   iflag1    9d12 + iflag1*10 +9
00107 *  ist +4   iflag2    9d12 + iflag2*10 +9
00108 *  ist +5   scamin    minimum y-scale
00109 *  ist +6   scamax    maximum y-scale
00110 *  ist +7   jdlast    address of the next histogram
00111 *                     from previous history of calls (see hadres)
00112 *          ----------------------------------
00113 *              Binning size informations
00114 *          ----------------------------------
00115 *  One dimensional histogram            Two dimensional histog.
00116 *  -------------------------            ----------------------
00117 *  (8-11) Binning information           (8-15) Binning information
00118 *  ist2 = ist+7
00119 *  ist2 +1    NCHX                          ist2 +5   NCHY
00120 *  ist2 +2      XL                          ist2 +6     YL
00121 *  ist2 +3      XU                          ist2 +7     YU
00122 *  ist2 +4   FACTX                          ist2 +8  FACTY
00123 *
00124 *          ----------------------------------
00125 *             All kind of sums and  maxwt
00126 *          ----------------------------------
00127 *  One dimensional histogram            Two dimensional histog.
00128 *  -------------------------            ----------------------
00129 *  (12-24) Under/over-flow average x    (16-24)
00130 *  ist3 = ist+11
00131 *  ist3 +1   Underflow                     All nine combinations
00132 *  ist3 +2   Normal                        (U,N,O) x (U,N,O)
00133 *  ist3 +3   Overflow                      sum wt only (no errors)
00134 *  ist3 +4   U  sum w**2
00135 *  ist3 +5   N  sum w**2
00136 *  ist3 +6   O  sum w**2
00137 *  ist3 +7   Sum 1
00138 *  ist3 +8   Sum wt*x
00139 *  ist3 +9   Sum wt*x*x
00140 *  ist3 +10  nevzer    (GLK_WtMon)
00141 *  ist3 +11  nevove    (GLK_WtMon)
00142 *  ist3 +12  nevacc    (GLK_WtMon)
00143 *  ist3 +13  maxwt     (GLK_WtMon)
00144 *          ----------------------------------
00145 *           Content of bins including errors
00146 *          ----------------------------------
00147 *  (25 to 24+2*nchx)                     (25 to 24 +nchx*nchy)
00148 *     sum wt and sum wt**2            sum wt only (no errors)
00149 *  ----------------------------------------------------------------
00150 *//////////////////////////////////////////////////////////////////////////////
00151 
00152       SUBROUTINE GLK_Initialize
00153 *     *************************
00154 * First Initialization called from may routines
00155 *     *************************************
00156       IMPLICIT NONE
00157 *----------------------------------------------------------------------
00158       INCLUDE 'GLK.h'
00159       SAVE
00160 *----------------------------------------------------------------------
00161 * Note that backslash definition is varying from one
00162 * instalation/compiler to another, you have to figure out by yourself
00163 * how to fill backslash code into m_BS
00164       CHARACTER*1 BBS1, BBS2
00165       DATA BBS1 /'\\'/    ! IBM or HP with 'f77 +B '
00166       DATA BBS2 /1H\ /    ! HP  f77 with standard options
00167 *-----------------------------------------------
00168       INTEGER init,i,k
00169       DATA init /0/
00170 *-----------------------------------------------
00171       IF(init .NE. 0) RETURN
00172       init=1
00173 * default output unit
00174       m_out=16
00175       m_length=0
00176 * color
00177       m_KeyCol=0
00178 * table range
00179       m_KeyTbr=0
00180       DO k=1,3
00181          m_TabRan(k)=1
00182       ENDDO
00183 *
00184       DO k=1,80
00185          m_Color(k:k)=' '
00186       ENDDO
00187       m_Color(1:1)='%'
00188 *
00189       DO i=1,m_idmax
00190          DO k=1,3
00191             m_index(i,k)=0
00192          ENDDO
00193          DO k=1,80
00194             m_titlc(i)(k:k)=' '
00195          ENDDO
00196       ENDDO
00197       DO k=1,m_LenmB
00198          m_b(k)=0d0
00199       ENDDO
00200 *----------------------------------------------------------------------
00201       m_BS = BBS1      ! IBM or HP with 'f77 +B '
00202 **    m_BS = BBS2      ! HP standard options
00203 *----------------------------------------------------------------------
00204       END
00205 
00206       SUBROUTINE GLK_Flush
00207 *     ********************
00208 * FLUSH memory, all histos erased!
00209 *     *************************************
00210       IMPLICIT NONE
00211       INCLUDE 'GLK.h'
00212       INTEGER i,k
00213 *------------------------------------------------
00214       CALL GLK_Initialize
00215       m_length=0
00216       DO i=1,m_idmax
00217          DO k=1,3
00218             m_index(i,k)=0
00219          ENDDO
00220          DO k=1,80
00221             m_titlc(i)(k:k)=' '
00222          ENDDO
00223       ENDDO
00224       DO k=1,m_LenmB
00225          m_b(k)=0d0
00226       ENDDO
00227       END
00228 
00229       LOGICAL FUNCTION GLK_Exist(id)
00230 *     ******************************
00231 * this function is true when id  exists !!!!
00232 *     ***************************
00233       IMPLICIT NONE
00234       INCLUDE 'GLK.h'
00235       INTEGER id,lact
00236 *------------------------------------------------
00237       CALL GLK_hadres(id,lact)
00238       GLK_Exist = lact .NE. 0
00239 ***   IF(GLK_Exist)      WRITE(6,*) 'GLK_Exist: does   ID,lact=',id,lact
00240 ***   IF(.not.GLK_Exist) write(6,*) 'GLK_Exist: doesnt ID,lact=',id,lact
00241       END
00242 
00243       DOUBLE PRECISION FUNCTION GLK_hi(id,ib)
00244 *     **********************
00245 * getting out bin content
00246 * S.J. 18-Nov. 90
00247 *     ***********************************
00248       IMPLICIT NONE
00249       INCLUDE 'GLK.h'
00250       INTEGER  id,ib
00251 * locals
00252       INTEGER    ist,ist2,ist3,iflag2,ityphi,nch,idmem,lact
00253       SAVE idmem
00254       DATA idmem / -1256765/
00255 *------------------------------------------------------
00256       IF(id .EQ. idmem) goto 100
00257       idmem=id
00258 * some checks, not repeated if id the same as previously
00259       CALL GLK_hadres(id,lact)
00260       IF(lact .EQ. 0) THEN
00261         CALL GLK_Stop1(' GLK_hi: nonexisting histo id=',id)
00262       ENDIF
00263       ist  = m_index(lact,2)
00264       ist2 = ist+7
00265       ist3 = ist+11
00266 * checking if histo is of proper type
00267       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00268       ityphi   = mod(iflag2,10)
00269       IF(ityphi .NE. 1 .AND. ityphi.NE.3) THEN
00270          CALL GLK_Stop1(' GLK_hi: 1-dim histos only !!! id=',id)
00271       ENDIF
00272   100 continue
00273       nch  = nint(m_b(ist2+1))
00274       IF(ib .EQ. 0) THEN
00275 * underflow
00276          GLK_hi=   m_b(ist3 +1)
00277       ELSEIF(ib .GE. 1.and.ib .LE. nch) THEN
00278 * normal bin
00279          GLK_hi=   m_b(ist +m_buf1+ib)
00280       ELSEIF(ib .EQ. nch+1) THEN
00281 * overflow
00282          GLK_hi=   m_b(ist3 +3)
00283       ELSE
00284 * abnormal exit
00285          CALL GLK_Stop1(' GLK_hi: wrong binning id=',id)
00286       ENDIF
00287       END
00288 
00289       DOUBLE PRECISION FUNCTION  GLK_hie(id,ib)
00290 *     ************************
00291 * getting out error of the bin
00292 * s.j. 18-nov. 90
00293 *     ***********************************
00294       IMPLICIT NONE
00295       INCLUDE 'GLK.h'
00296 * locals
00297       INTEGER    ist,ist2,ist3,iflag2,ityphi,nch,lact,ib,id
00298       SAVE       idmem
00299       INTEGER    idmem
00300       DATA idmem / -1256765/
00301 *---------------------------------------------------------
00302       IF(id .EQ. idmem) goto 100
00303       idmem=id
00304 * some checks, not repeated if id the same as previously
00305       CALL GLK_hadres(id,lact)
00306       IF(lact .EQ. 0) THEN
00307         CALL GLK_Stop1(' GLK_hie: nonexisting histo id=',id)
00308       ENDIF
00309       ist  = m_index(lact,2)
00310       ist2 = ist+7
00311       ist3 = ist+11
00312 * checking if histo is of proper type
00313       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00314       ityphi   = mod(iflag2,10)
00315       IF(ityphi .NE. 1) THEN
00316         CALL GLK_Stop1(' GLK_hie: 1-dim histos only !!! id=',id)
00317       ENDIF
00318   100 CONTINUE
00319       nch  = m_b(ist2+1)
00320       IF(ib .EQ. 0) THEN
00321 * underflow
00322          GLK_hie=   dsqrt( dabs(m_b(ist3 +4)))
00323       ELSEIF(ib .GE. 1.and.ib .LE. nch) THEN
00324 *...normal bin, error content
00325          GLK_hie=   dsqrt( dabs(m_b(ist+m_buf1+nch+ib)) )
00326       ELSEIF(ib .EQ. nch+1) THEN
00327 * overflow
00328          GLK_hie=   dsqrt( dabs(m_b(ist3 +6)))
00329       ELSE
00330 * abnormal exit
00331          CALL GLK_Stop1('+++GLK_hie: wrong binning id= ',id)
00332       ENDIF
00333       END
00334 
00335       SUBROUTINE GLK_Fil1(id,xx,wtx)
00336 *     *****************************
00337 * recommended fast filling 1-dim. histogram s.j. 18 nov. 90
00338 * overflow/underflow corrected by Maciek and Zbyszek
00339 *     ***********************************
00340       IMPLICIT NONE
00341       INCLUDE 'GLK.h'
00342       INTEGER           id
00343       DOUBLE PRECISION  xx,wtx
00344 * local
00345       INTEGER           ist,ist2,ist3,iflag2,ityphi,ipose1,iposx1,kposx1,kpose1,kx,nchx,lact
00346       DOUBLE PRECISION  x1,wt1,xl,factx,xu
00347 *-------------------------------------------------------------------------
00348       CALL GLK_hadres(id,lact)
00349 * exit for non-existig histo
00350       IF(lact .EQ. 0)  RETURN
00351       ist  = m_index(lact,2)
00352       ist2 = ist+7
00353       ist3 = ist+11
00354 * one-dim. histo only
00355       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00356       ityphi   = mod(iflag2,10)
00357       IF(ityphi .NE. 1)  CALL GLK_Stop1('+++GLK_Fil1: wrong id=  ',id)
00358       x1= xx
00359       wt1= wtx
00360       m_index(lact,3)=m_index(lact,3)+1
00361 * all entries
00362       m_b(ist3 +7)  =m_b(ist3 +7)   +1
00363 * for average x
00364       m_b(ist3 +8)  =m_b(ist3 +8)  +wt1
00365       m_b(ist3 +9)  =m_b(ist3 +9)  +wt1*x1
00366 * filling coordinates
00367       nchx  =m_b(ist2 +1)
00368       xl    =m_b(ist2 +2)
00369       xu    =m_b(ist2 +3)
00370       factx =m_b(ist2 +4)
00371       IF(x1 .LT. xl) THEN
00372 * underflow
00373          iposx1 = ist3 +1
00374          ipose1 = ist3 +4
00375          kposx1 = 0
00376       ELSEIF(x1 .GT. xu) THEN
00377 * or overflow
00378          iposx1 = ist3 +3
00379          ipose1 = ist3 +6
00380          kposx1 = 0
00381       ELSE
00382 * or any normal bin
00383          iposx1 = ist3 +2
00384          ipose1 = ist3 +5
00385 * or given normal bin
00386          kx = (x1-xl)*factx+1d0
00387          kx = MIN( MAX(kx,1) ,nchx)
00388          kposx1 = ist +m_buf1+kx
00389          kpose1 = ist +m_buf1+nchx+kx
00390       ENDIF
00391       m_b(iposx1) = m_b(iposx1)  +wt1
00392       m_b(ipose1) = m_b(ipose1)  +wt1*wt1
00393       IF( kposx1 .NE. 0) m_b(kposx1) = m_b(kposx1)  +wt1
00394       IF( kposx1 .NE. 0) m_b(kpose1) = m_b(kpose1)  +wt1*wt1
00395       END   !GLK_Fil1
00396 
00397       SUBROUTINE GLK_Fil1diff(id,xx,wtx,yy,wty)
00398 *     *****************************************
00399 * Special filling routine to fill the difference f(x)-g(y)
00400 * in the case when f and g are very similar x and y are close for each event.
00401 * In this case coherent filling is done if x and y fall into the same bin.
00402 * Note that bin width starts to be important in this method.
00403 *     ***********************************
00404       IMPLICIT NONE
00405       INCLUDE 'GLK.h'
00406 *
00407       INTEGER           id
00408       DOUBLE PRECISION  xx,wtx,yy,wty
00409 *    
00410       DOUBLE PRECISION  x1,x2,wt2,wt1,factx,xl,xu
00411       INTEGER           ist,ist2,ist3,iflag2,ityphi,kx,ke1,ie1,kx1,kx2,ke2,ix2,ie2,nchx,lact,ix1
00412 *-----------------------------------------------------------------
00413       CALL GLK_hadres(id,lact)
00414 * exit for non-existig histo
00415       IF(lact .EQ. 0)  RETURN
00416       ist  = m_index(lact,2)
00417       ist2 = ist+7
00418       ist3 = ist+11
00419 * one-dim. histo only
00420       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00421       ityphi   = mod(iflag2,10)
00422       IF(ityphi .NE. 1) THEN
00423         CALL GLK_Stop1('GLK_Fil1diff: 1-dim histos only !!! id=',id)
00424       ENDIF
00425       x1= xx
00426       x2= yy
00427       wt1= wtx
00428       wt2= wty
00429       m_index(lact,3)=m_index(lact,3)+1
00430 * all entries
00431       m_b(ist3 +7)  =m_b(ist3 +7)   +1
00432 * for average x or y not very well defined yet
00433       m_b(ist3 +8)  =m_b(ist3 +8)  +wt1*x1 - wt2*x2
00434       m_b(ist3 +9)  =m_b(ist3 +9)  +wt1*x1*x1 - wt2*x2*x2
00435 * filling coordinates
00436       nchx  =m_b(ist2 +1)
00437       xl    =m_b(ist2 +2)
00438       xu    =m_b(ist2 +3)
00439       factx =m_b(ist2 +4)
00440 * first variable
00441       IF(x1 .LT. xl) THEN       ! underflow
00442          ix1 = ist3 +1
00443          ie1 = ist3 +4
00444          kx1 = 0
00445       ELSEIF(x1 .GT. xu) THEN   ! or overflow
00446          ix1 = ist3 +3
00447          ie1 = ist3 +6
00448          kx1 = 0
00449       ELSE                      ! normal bin
00450          ix1 = ist3 +2
00451          ie1 = ist3 +5
00452          kx = (x1-xl)*factx+1d0
00453          kx = MIN( MAX(kx,1) ,nchx)
00454          kx1 = ist +m_buf1+kx
00455          ke1 = ist +m_buf1+nchx+kx
00456       ENDIF
00457 * second variable
00458       IF(x2 .LT. xl) THEN       ! underflow
00459          ix2 = ist3 +1
00460          ie2 = ist3 +4
00461          kx2 = 0
00462       ELSEIF(x2 .GT. xu) THEN   ! or overflow
00463          ix2 = ist3 +3
00464          ie2 = ist3 +6
00465          kx2 = 0
00466       ELSE                      ! normal bin
00467          ix2 = ist3 +2
00468          ie2 = ist3 +5
00469          kx = (x2-xl)*factx+1d0
00470          kx = MIN( MAX(kx,1) ,nchx)
00471          kx2 = ist +m_buf1+kx
00472          ke2 = ist +m_buf1+nchx+kx
00473       ENDIF
00474 * coherent filling
00475       IF( ix1 .EQ. ix2 ) THEN
00476          m_b(ix1) = m_b(ix1)  +wt1-wt2
00477          m_b(ie1) = m_b(ie1)  +(wt1-wt2)**2
00478       ELSE
00479          m_b(ix1) = m_b(ix1)  +wt1
00480          m_b(ie1) = m_b(ie1)  +wt1*wt1
00481          m_b(ix2) = m_b(ix2)  -wt2
00482          m_b(ie2) = m_b(ie2)  +wt2*wt2
00483       ENDIF
00484       IF( kx1 .EQ. kx2 ) THEN
00485          IF( kx1 .NE. 0) THEN
00486             m_b(kx1) = m_b(kx1)  +wt1-wt2
00487             m_b(ke1) = m_b(ke1)  +(wt1-wt2)**2
00488          ENDIF
00489       ELSE
00490          IF( kx1 .NE. 0) THEN
00491             m_b(kx1) = m_b(kx1)  +wt1
00492             m_b(ke1) = m_b(ke1)  +wt1*wt1
00493          ENDIF
00494          IF( kx2 .NE. 0) THEN
00495             m_b(kx2) = m_b(kx2)  -wt2
00496             m_b(ke2) = m_b(ke2)  +wt2*wt2
00497          ENDIF
00498       ENDIF
00499       END   !GLK_Fil1diff
00500 
00501       SUBROUTINE GLK_Fil2(id,x,y,wtw)
00502 *     ****************************
00503 * this routine not finished, 1-dim only!
00504 *     ***********************************
00505       IMPLICIT NONE     
00506       INCLUDE 'GLK.h'
00507       INTEGER           id
00508       DOUBLE PRECISION  x,y,wtw
00509 * local
00510       INTEGER           ist,iflag2,ityphi,ist2,ist3,nchx,nchy,ly,ky,k2,kx,lact,lx,k,l
00511       DOUBLE PRECISION  xx,yy,wt,factx,xl,yl,facty
00512 *-------------------------------------------------------
00513       CALL GLK_hadres(id,lact)
00514       IF(lact .EQ. 0)  RETURN
00515       ist  = m_index(lact,2)
00516 * one-dim. histo
00517       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00518       ityphi   = mod(iflag2,10)
00519       IF(ityphi .NE. 2) THEN
00520         CALL GLK_Stop1('GLK_Fil2: 2-dim histos only !!! id=',id)
00521       ENDIF
00522 *...two-dim. scattergram, no errors!
00523       ist2 = ist+7
00524       ist3 = ist+15
00525       xx= x
00526       yy= y
00527       wt= wtw
00528       m_index(lact,3)=m_index(lact,3)+1
00529 * x-axis
00530       nchx  =m_b(ist2 +1)
00531       xl    =m_b(ist2 +2)
00532       factx =m_b(ist2 +4)
00533       kx=(xx-xl)*factx+1d0
00534       lx=2
00535       IF(kx .LT. 1)     lx=1
00536       IF(kx .GT. nchx)  lx=3
00537       l     = ist+34  +lx
00538       m_b(l)  = m_b(l)    +wt
00539       k     = ist+m_buf2  +kx
00540       IF(lx .EQ. 2) m_b(k)  =m_b(k)  +wt
00541       k2    = ist+m_buf2  +nchx+kx
00542       IF(lx .EQ. 2) m_b(k2) =m_b(k2) +wt**2
00543 * y-axix
00544       nchy  =m_b(ist2 +5)
00545       yl    =m_b(ist2 +6)
00546       facty =m_b(ist2 +8)
00547       ky=(yy-yl)*facty+1d0
00548       ly=2
00549       IF(ky .LT. 1)    ly=1
00550       IF(ky .GT. nchy) ly=3
00551 * under/over-flow
00552       l = ist3  +lx +3*(ly-1)
00553       m_b(l) =m_b(l)+wt
00554 * regular bin
00555       k = ist+m_buf2 +kx +nchx*(ky-1)
00556       IF(lx .EQ. 2.and.ly .EQ. 2) m_b(k)=m_b(k)+wt
00557       END
00558 
00559       SUBROUTINE GLK_Book1(id,title,nnchx,xxl,xxu)
00560 *     ********************************************
00561       IMPLICIT NONE
00562       INCLUDE 'GLK.h'
00563       INTEGER          id
00564       DOUBLE PRECISION xxl,xxu
00565       CHARACTER*80 title
00566 * locals
00567       DOUBLE PRECISION xl,xu,ddx
00568       INTEGER          ist,nchx,ioplog,iopsla,ioperb,iflag2,ityphi,iflag1
00569       INTEGER          ist3,ist2,lengt2,lact,nnchx,iopsc2,iopsc1,j
00570       LOGICAL GLK_Exist
00571 *-------------------------------------------------
00572       CALL GLK_Initialize
00573       IF(GLK_Exist(id)) goto 900
00574       ist=m_length
00575       CALL GLK_hadres(0,lact)
00576 * Check if there is a free entry in the m_index
00577       IF(lact .EQ. 0)
00578      $     CALL GLK_Stop1('GLK_Book1: to many histos !!!!!,   id= ',id)
00579       m_index(lact,1)=id
00580       m_index(lact,2)=m_length
00581       m_index(lact,3)=0
00582 * -------
00583       CALL GLK_Copch(title,m_titlc(lact))
00584       nchx =nnchx
00585       IF(nchx .GT. m_MaxNb)
00586      $     CALL GLK_Stop1(' GLK_Book1: To many bins requested,id= ',id)
00587       xl   =xxl
00588       xu   =xxu
00589 * ---------- title and bin content ----------
00590       lengt2 = m_length +2*nchx +m_buf1+1
00591       IF(lengt2 .GE. m_LenmB)
00592      $  CALL GLK_Stop1('GLK_Book1:too litle storage, m_LenmB= ',m_LenmB)
00593 *
00594       DO j=m_length+1,lengt2+1
00595          m_b(j) = 0d0
00596       ENDDO
00597       m_length=lengt2
00598 *... default flags
00599       ioplog   = 1
00600       iopsla   = 1
00601       ioperb   = 1
00602       iopsc1   = 1
00603       iopsc2   = 1
00604       iflag1   =
00605      $ ioplog+10*iopsla+100*ioperb+1000*iopsc1+10000*iopsc2
00606       ityphi   = 1
00607       iflag2   = ityphi
00608 * examples of decoding flags
00609 *      id       = nint(m_b(ist+2)-9d0-9d12)/10
00610 *      iflag1   = nint(m_b(ist+3)-9d0-9d12)/10
00611 *      ioplog = mod(iflag1,10)
00612 *      iopsla = mod(iflag1,100)/10
00613 *      ioperb = mod(iflag1,1000)/100
00614 *      iopsc1 = mod(iflag1,10000)/1000
00615 *      iopsc2 = mod(iflag1,100000)/10000
00616 *      iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
00617 *      ityphi = mod(iflag2,10)
00618 *--------- buffer -----------------
00619 * header
00620       m_b(ist +1)  = 9999999999999d0
00621       m_b(ist +2)  = 9d12 +     id*10 +9d0
00622       m_b(ist +3)  = 9d12 + iflag1*10 +9d0
00623       m_b(ist +4)  = 9d12 + iflag2*10 +9d0
00624 * dummy vertical scale
00625       m_b(ist +5)  =  -100d0
00626       m_b(ist +6)  =   100d0
00627 * pointer used to speed up search of histogram address
00628       m_b(ist +7)  =   0d0
00629 * information on binning
00630       ist2         = ist+7
00631       m_b(ist2 +1) = nchx
00632       m_b(ist2 +2) = xl
00633       m_b(ist2 +3) = xu
00634       ddx = xu-xl
00635       IF(ddx .EQ. 0d0) CALL GLK_Stop1('+++GLK_Book1: xl=xu,  id= ',id)
00636       m_b(ist2 +4) = DFLOAT(nchx)/ddx
00637 *
00638 * under/over-flow etc.
00639       ist3       = ist+11
00640       DO j=1,13
00641          m_b(ist3 +j)=0d0
00642       ENDDO
00643       RETURN
00644 *----------------
00645  900  CALL GLK_Retu1(' WARNING GLK_Book1: already exists id= ', id)
00646       END
00647 
00648       SUBROUTINE GLK_Retu1(mesage,id)
00649 *     *******************************
00650       IMPLICIT NONE
00651       INCLUDE 'GLK.h'
00652       SAVE
00653       INTEGER id
00654       CHARACTER*(*) mesage
00655 *-----------------------------
00656       WRITE(m_out,'(a)')
00657      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00658       WRITE(m_out,'(a,a,i10,a)')
00659      $                          '++ ', mesage, id, ' ++'
00660       WRITE(m_out,'(a)')
00661      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00662       WRITE(6   ,'(a)')
00663      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00664       WRITE(6   ,'(a,a,i10,a)')
00665      $                          '++ ', mesage, id, ' ++'
00666       WRITE(6   ,'(a)')
00667      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00668       END
00669 
00670       SUBROUTINE GLK_Stop1(mesage,id)
00671 *     *******************************
00672       IMPLICIT NONE
00673       INCLUDE 'GLK.h'
00674       SAVE
00675       CHARACTER*(*) mesage
00676       INTEGER id
00677 *-----------------------------
00678       WRITE(m_out,'(a)')
00679      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00680       WRITE(m_out,'(a,a,i10,a)')
00681      $                          '++ ', mesage, id, ' ++'
00682       WRITE(m_out,'(a)')
00683      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00684       WRITE(6   ,'(a)')
00685      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00686       WRITE(6   ,'(a,a,i10,a)')
00687      $                          '++ ', mesage, id, ' ++'
00688       WRITE(6   ,'(a)')
00689      $          '++++++++++++++++++++++++++++++++++++++++++++++++++++++'
00690       STOP
00691       END
00692 
00693 
00694       SUBROUTINE GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
00695 *     ********************************************************
00696 * decoding option flags
00697 *     **********************************
00698       IMPLICIT NONE
00699       INCLUDE 'GLK.h'
00700       INTEGER    id,ioplog,iopsla,ioperb,iopsc1,iopsc2
00701       INTEGER    ist,iflag1,lact
00702 *----------------------------------------------------------------
00703       CALL GLK_hadres(id,lact)
00704       IF(lact .EQ. 0) RETURN
00705       ist=m_index(lact,2)
00706 * decoding flags
00707       iflag1   = nint(m_b(ist+3)-9d0-9d12)/10
00708       ioplog = mod(iflag1,10)
00709       iopsla = mod(iflag1,100)/10
00710       ioperb = mod(iflag1,1000)/100
00711       iopsc1 = mod(iflag1,10000)/1000
00712       iopsc2 = mod(iflag1,100000)/10000
00713       END
00714 
00715       SUBROUTINE GLK_idopt(id,ch)
00716 *     ************************
00717       IMPLICIT NONE
00718       INCLUDE 'GLK.h'
00719       INTEGER     id
00720       CHARACTER*4 ch
00721 *
00722       INTEGER     lact,ist,ioplog,ioperb,iopsla,iopsc1,iopsc2,iflag1
00723 *----------------------------------------------------------------
00724       CALL GLK_hadres(id,lact)
00725       IF(lact .EQ. 0) RETURN
00726       ist=m_index(lact,2)
00727 * decoding flags
00728       CALL GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
00729       IF(ch .EQ.       'LOGY'  ) THEN
00730 * log scale for print
00731         ioplog = 2
00732       ELSEIF(ch .EQ.   'ERRO'  ) THEN
00733 * errors in printing/plotting
00734        ioperb  = 2
00735       ELSEIF(ch .EQ.   'SLAN'  ) THEN
00736 * slanted line in plotting
00737        iopsla  = 2
00738       ELSEIF(ch .EQ.   'YMIN'  ) THEN
00739        iopsc1  = 2
00740       ELSEIF(ch .EQ.   'YMAX'  ) THEN
00741        iopsc2  = 2
00742       ENDIF
00743 * encoding back
00744       iflag1   = ioplog+10*iopsla+100*ioperb+1000*iopsc1+10000*iopsc2
00745       m_b(ist+3) = 9d12 + iflag1*10 +9d0
00746       END
00747 
00748 
00749       SUBROUTINE GLK_BookFun1(id,title,nchx,xmin,xmax,func)
00750 */////////////////////////////////////////////////////////////////////////
00751 *//   fills histogram with function func(x)                             //
00752 */////////////////////////////////////////////////////////////////////////
00753       IMPLICIT NONE
00754       INCLUDE 'GLK.h'
00755       INTEGER           id
00756       DOUBLE PRECISION  xmin,xmax,func
00757       CHARACTER*80 title
00758 *
00759       DOUBLE PRECISION  yy(m_MaxNb)
00760       EXTERNAL func
00761       LOGICAL GLK_Exist
00762       INTEGER           ib,nchx
00763       DOUBLE PRECISION  xl,xu,x
00764 *---------------------------------------------------------------------
00765       CALL GLK_Initialize
00766       IF(GLK_Exist(id)) GOTO 900
00767  15   xl=xmin
00768       xu=xmax
00769       CALL GLK_Book1(id,title,nchx,xl,xu)
00770 *...slanted line in plotting
00771       CALL GLK_idopt(id,'SLAN')
00772       IF(nchx .GT. 200) goto 901
00773       DO ib=1,nchx
00774          x= xmin +(xmax-xmin)/nchx*(ib-0.5d0)
00775          yy(ib) = func(x)
00776       ENDDO
00777       CALL GLK_Pak(id,yy)
00778       RETURN
00779  900  CALL GLK_Retu1('+++GLK_BookFun1: already exists id=',id)
00780       CALL GLK_Delet(id)
00781       GOTO 15
00782  901  CALL GLK_Stop1('+++GLK_BookFun1: to many bins, id=',id)
00783       END
00784 
00785       SUBROUTINE GLK_BookFun1I(id,title,nchx,xmin,xmax,func)
00786 */////////////////////////////////////////////////////////////////////////
00787 *//   Fills histogram with function func(x)                             //
00788 *//   Gauss integration over each bin is done, can be slow.             //
00789 */////////////////////////////////////////////////////////////////////////
00790       IMPLICIT NONE
00791       INCLUDE 'GLK.h'
00792       INTEGER           id
00793       DOUBLE PRECISION  xmin,xmax,func
00794       CHARACTER*80 title
00795 *
00796       DOUBLE PRECISION  yy(m_MaxNb)
00797       EXTERNAL func
00798       LOGICAL GLK_Exist
00799       INTEGER           ib,nchx
00800       DOUBLE PRECISION  xl,xu,x
00801       DOUBLE PRECISION  GLK_Gauss,a,b,Eeps,dx
00802 *---------------------------------------------------------------------
00803       CALL GLK_Initialize
00804       IF(GLK_Exist(id)) GOTO 900
00805  15   xl=xmin
00806       xu=xmax
00807       CALL GLK_Book1(id,title,nchx,xl,xu)
00808       IF(nchx .GT. 200) goto 901
00809       Eeps = -0.01d0             !!! relat. precision requirement not very demanding
00810       dx = (xmax-xmin)/nchx
00811       DO ib=1,nchx
00812          a= xmin +dx*(ib-1)
00813          b= xmin +dx*ib
00814          yy(ib) = GLK_Gauss(func,a,b,Eeps)/dx   !! 16-point Gauss integration over bin
00815       ENDDO
00816       CALL GLK_Pak(id,yy)
00817       RETURN
00818  900  CALL GLK_Retu1('+++GLK_BookFun1I: already exists id=',id)
00819       CALL GLK_Delet(id)
00820       GOTO 15
00821  901  CALL GLK_Stop1('+++GLK_BookFun1I: to many bins, id=',id)
00822       END
00823 
00824       SUBROUTINE GLK_BookFun1S(id,title,nchx,xmin,xmax,func)
00825 */////////////////////////////////////////////////////////////////////////
00826 *// Fills histogram with function func(x)                               //
00827 *// three point fit used                                                //
00828 */////////////////////////////////////////////////////////////////////////
00829       IMPLICIT NONE
00830       INCLUDE 'GLK.h'
00831       DOUBLE PRECISION  xmin,xmax,func
00832       EXTERNAL func
00833       INTEGER  id,nchx
00834       CHARACTER*80 title
00835 * locals
00836       DOUBLE PRECISION  yy(m_MaxNb),yy1(0:m_MaxNb)
00837       LOGICAL GLK_Exist
00838       DOUBLE PRECISION  xl,xu,x3,x2,dx
00839       INTEGER           ib  
00840 *--------------------------------------------------------
00841       CALL GLK_Initialize
00842       IF( GLK_Exist(id) ) GOTO 900
00843  15   xl=xmin
00844       xu=xmax
00845       CALL GLK_Book1(id,title,nchx,xl,xu)
00846 
00847 *...slanted line in plotting
00848       CALL GLK_idopt(id,'SLAN')
00849       IF(nchx.gt.200) GOTO 901
00850 
00851       yy1(0) = func(xmin)
00852       dx=(xmax-xmin)/nchx
00853 
00854       DO ib=1,nchx
00855          x2= xmin +dx*(ib-0.5d0)
00856          x3= x2 +dx*0.5d0
00857          yy(ib)  = func(x2)
00858          yy1(ib) = func(x3)
00859 *..  simpson
00860          yy(ib) = ( yy1(ib-1) +4*yy (ib) +yy1(ib))/6d0
00861       ENDDO
00862 
00863       CALL GLK_Pak(id,yy)
00864       RETURN
00865  900  CALL GLK_Retu1('+++GLK_BookFun1S: already exists, id=',id)
00866       CALL GLK_Delet(id)
00867       GOTO 15
00868  901  CALL GLK_Stop1(' +++GLK_BookFun1S: to many bins, id=',id)
00869       END
00870 
00871       DOUBLE PRECISION FUNCTION GLK_Gauss(f,a,b,Eeps)
00872 *//////////////////////////////////////////////////////////////////////////////
00873 *//                                                                          //
00874 *// This is iterative integration procedure                                  //
00875 *// originates  probably from CERN library                                   //
00876 *// it subdivides inegration range until required precision is reached       //
00877 *// precision is a difference from 8 and 16 point Gauss itegr. result        //
00878 *// Eeps POSITIVE treated as absolute precision                              //
00879 *// Eeps NEGATIVE treated as relative precision                              //
00880 *//                                                                          //
00881 *//////////////////////////////////////////////////////////////////////////////
00882       IMPLICIT NONE
00883       DOUBLE PRECISION  f,a,b,Eeps
00884 *
00885       DOUBLE PRECISION  c1,c2,bb,s8,s16,y,aa,const,delta,eps,u
00886       INTEGER           i
00887 *
00888       DOUBLE PRECISION  w(12),x(12)
00889       EXTERNAL f
00890       DATA const /1.0d-19/
00891       DATA w
00892      1/0.10122 85362 90376, 0.22238 10344 53374, 0.31370 66458 77887,
00893      2 0.36268 37833 78362, 0.02715 24594 11754, 0.06225 35239 38648,
00894      3 0.09515 85116 82493, 0.12462 89712 55534, 0.14959 59888 16577,
00895      4 0.16915 65193 95003, 0.18260 34150 44924, 0.18945 06104 55069/
00896       DATA x
00897      1/0.96028 98564 97536, 0.79666 64774 13627, 0.52553 24099 16329,
00898      2 0.18343 46424 95650, 0.98940 09349 91650, 0.94457 50230 73233,
00899      3 0.86563 12023 87832, 0.75540 44083 55003, 0.61787 62444 02644,
00900      4 0.45801 67776 57227, 0.28160 35507 79259, 0.09501 25098 37637/
00901 *-----------------------------------------------------------------------------
00902       eps=abs(Eeps)
00903       delta=const*abs(a-b)
00904       GLK_Gauss=0d0
00905       aa=a
00906     5 y=b-aa
00907       IF(abs(y)  .LE.  delta) RETURN
00908     2 bb=aa+y
00909       c1=0.5d0*(aa+bb)
00910       c2=c1-aa
00911       s8=0d0
00912       s16=0d0
00913       DO 1 i=1,4
00914       u=x(i)*c2
00915     1 s8=s8+w(i)*(f(c1+u)+f(c1-u))
00916       DO 3 i=5,12
00917       u=x(i)*c2
00918     3 s16=s16+w(i)*(f(c1+u)+f(c1-u))
00919       s8=s8*c2
00920       s16=s16*c2
00921       IF(Eeps .LT. 0d0) THEN
00922         IF(abs(s16-s8)  .GT.  eps*abs(s16)) GOTO 4
00923       ELSE
00924         IF(abs(s16-s8)  .GT.  eps) GOTO 4
00925       ENDIF
00926       GLK_Gauss=GLK_Gauss+s16
00927       aa=bb
00928       GOTO 5
00929     4 y=0.5d0*y
00930       IF(abs(y)  .GT.  delta) GOTO 2
00931       WRITE(*,7)
00932       GLK_Gauss=0d0
00933       RETURN
00934     7 FORMAT(1x,36hgaus  ... too high accuracy required)
00935       END
00936 
00937 
00938 
00939       SUBROUTINE GLK_Book2(ID,TITLE,NCHX,XL,XU,NCHY,YL,YU)
00940 *     ***************************************************
00941       IMPLICIT NONE
00942       INCLUDE 'GLK.h'
00943       INTEGER           ID,NCHX,NCHY
00944       DOUBLE PRECISION  XL,XU,YL,YU
00945       CHARACTER*80 TITLE
00946 *
00947       INTEGER   ist,lact,lengt2,j,nnchx,nnchy
00948       LOGICAL GLK_EXIST
00949 *-------------------------------------------------------------------------
00950       CALL GLK_Initialize
00951       IF(GLK_EXIST(ID)) GOTO 900
00952       ist=m_length
00953       CALL GLK_hadres(0,lact)
00954       IF(LACT .EQ. 0) GOTO 901
00955       m_index(LACT,1)=ID
00956       m_index(LACT,2)=m_length
00957       CALL GLK_COPCH(TITLE,M_TITLC(LACT))
00958       nnchx=NCHX
00959       nnchy=NCHY
00960       LENGT2 = M_LENGTH  +44+nnchx*nnchy
00961       IF(LENGT2 .GE. m_LenmB) GOTO 902
00962       DO 10 J=M_LENGTH+1,LENGT2+1
00963    10 m_b(J) = 0D0
00964       M_LENGTH=LENGT2
00965       m_b(ist+1)=nnchx
00966       m_b(ist+2)=XL
00967       m_b(ist+3)=XU
00968       m_b(ist+4)=float(nnchx)/(m_b(ist+3)-m_b(ist+2))
00969       m_b(ist+5)=nnchy
00970       m_b(ist+6)=YL
00971       m_b(ist+7)=YU
00972       m_b(ist+8)=float(nnchy)/(m_b(ist+7)-m_b(ist+6))
00973       RETURN
00974   900 CALL GLK_Retu1('GLK_Book2: histo already exists!!!! id=',id)
00975       RETURN
00976   901 CALL GLK_Stop1('GLK_Book2: too many histos !!!!! lact= ',LACT)
00977       RETURN
00978   902 CALL GLK_Stop1('GLK_Book2: too litle storage, m_LenmB=',m_LenmB)
00979       RETURN
00980       END
00981 
00982       SUBROUTINE GLK_PrintAll
00983 *     ***********************
00984       IMPLICIT NONE
00985       INCLUDE 'GLK.h'
00986       SAVE
00987       INTEGER i,id
00988 
00989       DO i=1,m_idmax
00990          id=m_index(i,1)
00991          IF(id .GT. 0) CALL GLK_Print(id)
00992       ENDDO
00993       END
00994 
00995       SUBROUTINE GLK_ListPrint(mout)
00996 *//////////////////////////////////////////////////////////////////////////////
00997 *//                                                                          //
00998 *//////////////////////////////////////////////////////////////////////////////
00999       IMPLICIT NONE
01000       INCLUDE 'GLK.h'
01001       INTEGER i,id
01002       CHARACTER*80 title
01003       INTEGER             nb,mout
01004       DOUBLE PRECISION    xmin,xmax
01005 *----------------------------------
01006       WRITE(mout,*) 
01007      $'============================================================================================'
01008       WRITE(mout,*) 
01009      $'     ID                 TITLE                                                nb, xmin, xmax'
01010       DO i=1,m_idmax
01011          id=m_index(i,1)
01012          IF(id .NE. 0) THEN
01013             CALL GLK_hinbo1(id,title,nb,xmin,xmax)
01014             WRITE(mout,'(i8,a,a,i8,2g14.6)') id, '  ', title, nb,xmin,xmax
01015          ENDIF
01016       ENDDO
01017       END
01018 
01019 
01020 
01021       SUBROUTINE GLK_Print(id)
01022 *     ***********************
01023       IMPLICIT NONE
01024       INCLUDE 'GLK.h'
01025       INTEGER  id
01026 *
01027       DOUBLE PRECISION xl,bind,xlow,z,er,avex,dx,fact,ovef,undf,bmax,bmin,deltb
01028       DOUBLE PRECISION sum,sumw,sumx
01029       INTEGER          ist,ist2,ist3,idec,k2,k1,kros,j,ind,i,n,i1,ky,nchy,kx,nent,iflag2,lmx
01030       INTEGER          ioplog,iopsla,ioperb,iopsc1,iopsc2,lact,ker,ityphi,kzer,k,ibn,nchx,istr
01031       LOGICAL llg
01032       CHARACTER*1 line(0:105),lchr(22),lb,lx,li,l0
01033       SAVE lb,lx,li,l0,lchr
01034       DATA lb,lx,li,l0 /' ','X','I','0'/
01035       DATA lchr/' ','1','2','3','4','5','6','7','8','9',
01036      $      'A','B','C','D','E','F','G','H','I','J','K','*'/
01037 *---------------------------------------------------------------------------------
01038       CALL GLK_hadres(id,lact)
01039       IF(lact .EQ. 0) goto 900
01040       ist  = m_index(lact,2)
01041       ist2 = ist+7
01042       ist3 = ist+11
01043       idec    = nint(m_b(ist+2)-9d0-9d12)/10
01044       IF(idec .NE. id) WRITE(6,*) '++++GLK_PRINT: PANIC! ID,IDEC= ',ID,IDEC
01045       CALL GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
01046       ker    =  ioperb-1
01047       lmx = 67
01048       IF(ker .EQ. 1) lmx=54
01049       nent=m_index(lact,3)
01050       IF(nent  .EQ.  0)                          GOTO 901
01051       WRITE(m_out,1000) id,m_titlc(lact)
01052  1000 FORMAT('1',/,1X,I9,10X,A)
01053 *
01054 * one-dim. histo
01055       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01056       ityphi   = mod(iflag2,10)
01057       IF(ityphi .EQ. 2) GOTO 200
01058       IF( (ityphi.NE.1) .AND. (ityphi.NE.3) )
01059      $   CALL GLK_Stop1(' GLK_PRINT: wrong histo type, id=',id)
01060 
01061       nchx =   m_b(ist2 +1)
01062       xl   =   m_b(ist2 +2)
01063       dx   =  (  m_b(ist2 +3)-m_b(ist2 +2)  )/float(nchx)
01064 * fixing vertical scale
01065       istr=ist+m_buf1+1
01066       bmin = m_b(istr)
01067       bmax = m_b(istr)+1d-5*abs(m_b(istr))  ! problems for single bin case
01068       DO  ibn=istr,istr+nchx-1
01069          bmax = max(bmax,m_b(ibn))
01070          bmin = min(bmin,m_b(ibn))
01071       ENDDO
01072       IF(bmin  .EQ.  bmax)                       GOTO 903
01073       IF(iopsc1 .EQ. 2) bmin=m_b(ist +5)
01074       IF(iopsc2 .EQ. 2) bmax=m_b(ist +6)
01075 *
01076       llg=ioplog .EQ. 2
01077       IF(llg.and.bmin .LE. 0d0) bmin=bmax/10000.d0
01078 *
01079       deltb = bmax-bmin
01080       IF(deltb  .EQ.  0d0)                       GOTO 902
01081       fact  = (lmx-1)/deltb
01082       kzer  = -bmin*fact+1.00001d0
01083       IF(llg) fact=(lmx-1)/(log(bmax)-log(bmin))
01084       IF(llg) kzer=-log(bmin)*fact+1.00001d0
01085 *
01086       undf = m_b(ist3 +1)
01087       ovef = m_b(ist3 +3)
01088       avex = 0d0
01089       sumw  = m_b(ist3 +8)
01090       sumx  = m_b(ist3 +9)
01091       IF(sumw .NE. 0d0) avex = sumx/sumw
01092       WRITE(m_out,'(4a15      )')  'nent',' sum','bmin','bmax'
01093       WRITE(m_out,'(i15,3e15.5)')   nent,   sum,  bmin,  bmax
01094       WRITE(m_out,'(4a15  )')      'undf','ovef','sumw','avex'
01095       WRITE(m_out,'(4e15.5)')       undf,  ovef,  sumw,  avex
01096 *
01097       IF(llg) write(m_out,1105)
01098  1105 format(35x,17hlogarithmic scale)
01099 *
01100       kzer=max0(kzer,0)
01101       kzer=min0(kzer,lmx)
01102       xlow=xl
01103       do 100 k=1,nchx
01104 * first fill with blanks
01105       do  45 j=1,105
01106    45 line(j)  =lb
01107 * THEN fill upper and lower boundry
01108       line(1)  =li
01109       line(lmx)=li
01110       ind=istr+k-1
01111       bind=m_b(ind)
01112       bind= max(bind,bmin)
01113       bind= min(bind,bmax)
01114       kros=(bind-bmin)*fact+1.0001d0
01115       IF(llg) kros=log(bind/bmin)*fact+1.0001d0
01116       k2=max0(kros,kzer)
01117       k2=min0(lmx,max0(1,k2))
01118       k1=min0(kros,kzer)
01119       k1=min0(lmx,max0(1,k1))
01120       do 50 j=k1,k2
01121    50 line(j)=lx
01122       line(kzer)=l0
01123       z=m_b(ind)
01124       IF(ker .NE. 1) THEN
01125         WRITE(m_out,'(a, f7.4,  a, d14.6,  132a1)')
01126      $             ' ', xlow,' ',     z,' ',(line(i),i=1,lmx)
01127       ELSE
01128         er=dsqrt(dabs(m_b(ind+nchx)))
01129         WRITE(m_out,'(a,f7.4,  a,d14.6,  a,d14.6, 132a1 )')
01130      $             ' ',xlow,' ',    z,' ',   er,' ',(line(i),i=1,lmx)
01131       ENDIF
01132       xlow=xlow+dx
01133   100 continue
01134       RETURN
01135 *//////////////////////////////////////////////////////////////////////
01136 *// two dimensional requires complete restoration and tests          //
01137 *//////////////////////////////////////////////////////////////////////
01138   200 continue
01139       nchx=m_b(ist+1)
01140       nchy=m_b(ist+5)
01141       WRITE(m_out,2000) (lx,i=1,nchy)
01142  2000 format(1h ,10x,2hxx,100a1)
01143       do 300 kx=1,nchx
01144       do 250 ky=1,nchy
01145       k=ist +m_buf2 +kx+nchx*(ky-1)
01146       N=m_b(K)+1.99999D0
01147       n=max0(n,1)
01148       n=min0(n,22)
01149       IF(DABS(m_b(k)) .LT. 1D-20) n=1
01150       line(ky)=lchr(n)
01151   250 continue
01152       line(nchy+1)=lx
01153       i1=nchy+1
01154       WRITE(m_out,2100) (line(i),i=1,i1)
01155  2100 format(1h ,10x,1hx,100a1)
01156   300 continue
01157       WRITE(m_out,2000) (lx,i=1,nchy)
01158       RETURN
01159  900  CALL GLK_Retu1('GLK_PRINT: nonexisting histo',id)
01160       RETURN
01161  901  CALL GLK_Retu1('   GLK_PRINT: nent.eq.0',ID)
01162       RETURN
01163  902  CALL GLK_Retu1('   GLK_PRINT: wrong plotting limits, id=',id)
01164       RETURN
01165  903  CALL GLK_Retu1('   GLK_PRINT: bmin.eq.bmax, id=',id)
01166       END
01167 
01168       SUBROUTINE GLK_Operat(ida,chr,idb,idc,coef1,coef2)
01169 *     **********************************************
01170       IMPLICIT NONE
01171       INCLUDE 'GLK.h'
01172       INTEGER           ida,idb,idc
01173       DOUBLE PRECISION  coef1,coef2
01174       CHARACTER*80      title
01175       CHARACTER*1       chr
01176 *
01177       DOUBLE PRECISION  xl,xu
01178       INTEGER           ista,ista2,ista3,ncha,iflag2a,ityphia,lactb
01179       INTEGER           k,j,nchc,istc2,istc3,i1,j2,j3,j1,i2,i3,istc,istb2,istb3,nchb
01180       INTEGER           lacta,id,istb,nchx,iflag2b,ityphib,lactc
01181 *----------------------------------------------------------
01182       CALL GLK_hadres(ida,lacta)
01183       IF(lacta .EQ. 0) RETURN
01184       ista  = m_index(lacta,2)
01185       ista2 = ista+7
01186       ista3 = ista+11
01187       ncha  = m_b(ista2+1)
01188 * check for type
01189       iflag2a   = nint(m_b(ista+4)-9d0-9d12)/10
01190       ityphia   = mod(iflag2a,10)
01191       IF(ityphia .NE. 1) CALL GLK_Stop1('GLK_Operat: 1-dim histos only, id=',id)
01192 *------------------
01193       CALL GLK_hadres(idb,lactb)
01194       IF(lactb .EQ. 0) RETURN
01195       istb  = m_index(lactb,2)
01196       istb2 = istb+7
01197       istb3 = istb+11
01198       nchb  = m_b(istb2+1)
01199       IF(nchb .NE. ncha) goto 900
01200 * check for type
01201       iflag2b   = nint(m_b(istb+4)-9d0-9d12)/10
01202       ityphib   = mod(iflag2b,10)
01203       IF(ityphib .NE. 1)  CALL GLK_Stop1('GLK_Operat: 1-dim histos only, id=',id)
01204 *------------------
01205       CALL GLK_hadres(idc,lactc)
01206       IF(lactc .EQ. 0) THEN
01207 * ...if nonexistent, histo idc is here defined
01208         CALL GLK_hinbo1(ida,title,nchx,xl,xu)
01209         CALL GLK_Book1(idc,title,nchx,xl,xu)
01210         CALL GLK_hadres(idc,lactc)
01211         istc  = m_index(lactc,2)
01212 *...option copied from ida
01213         m_b(istc+ 3)= m_b(ista +3)
01214       ENDIF
01215 *...one nominal entry recorded
01216       m_index(lactc,3) = 1
01217 *
01218       istc  =  m_index(lactc,2)
01219       istc2 =  istc+7
01220       istc3 =  istc+11
01221       nchc  =  m_b(istc2+1)
01222 *
01223       IF(nchc .NE. ncha) goto 900
01224       IF(ncha .NE. nchb .OR. nchb .NE. nchc) GOTO 900
01225       DO k=1,ncha+3
01226          IF(k .GT. ncha) THEN
01227 *     underflow, overflow
01228             j=k-ncha
01229             i1 = ista3 +j
01230             i2 = istb3 +j
01231             i3 = istc3 +j
01232             j1 = ista3 +3+j
01233             j2 = istb3 +3+j
01234             j3 = istc3 +3+j
01235          ELSE
01236 *     normal bins
01237             i1 = ista +m_buf1 +k
01238             i2 = istb +m_buf1 +k
01239             i3 = istc +m_buf1 +k
01240             j1 = ista +m_buf1 +ncha+k
01241             j2 = istb +m_buf1 +ncha+k
01242             j3 = istc +m_buf1 +ncha+k
01243          ENDIF
01244          IF    (chr .EQ. '+')   THEN
01245             m_b(i3) =    coef1*m_b(i1) +    coef2*m_b(i2)
01246             m_b(j3) = coef1**2*m_b(j1) + coef2**2*m_b(j2)
01247          ELSEIF(chr .EQ. '-')   THEN
01248             m_b(i3) = coef1*m_b(i1) - coef2*m_b(i2)
01249             m_b(j3) = coef1**2*m_b(j1) + coef2**2*m_b(j2)
01250          ELSEIF(chr .EQ. '*')   THEN
01251             m_b(j3) = (coef1*coef2)**2
01252      $           *(m_b(j1)*m_b(i2)**2 + m_b(j2)*m_b(i1)**2)
01253             m_b(i3) = coef1*m_b(i1) * coef2*m_b(i2)
01254          ELSEIF(chr .EQ. '/')   THEN
01255             IF(m_b(i2) .EQ. 0d0) THEN
01256                m_b(i3) = 0d0
01257                m_b(j3) = 0d0
01258             ELSE
01259 ***               m_b(j3) = (coef1/coef2)**2/m_b(i2)**4           ! problems with overflow
01260 ***     $              *(m_b(j1)*m_b(i2)**2 + m_b(j2)*m_b(i1)**2) ! problems with overflow
01261                m_b(j3) = (coef1/coef2)**2 *m_b(j1) /m_b(i2)**2
01262      $                  +(coef1/coef2)**2 *m_b(j2) *(m_b(i1)/m_b(i2)**2)**2
01263                m_b(i3) = (coef1*m_b(i1) )/( coef2*m_b(i2))
01264             ENDIF
01265          ELSE
01266             GOTO 901
01267          ENDIF
01268       ENDDO
01269       RETURN
01270  900  WRITE(m_out,*) '+++++ GLK_Operat: non-equal no. bins ',ida,idb,idc
01271       WRITE(    6,*) '+++++ GLK_Operat: non-equal no. bins ',ida,idb,idc
01272       STOP
01273  901  WRITE(m_out,*) '+++++ GLK_Operat: wrong chr=',chr
01274       WRITE(    6,*) '+++++ GLK_Operat: wrong chr=',chr
01275       STOP
01276       END
01277 
01278       SUBROUTINE GLK_hinbo1(id,title,nchx,xl,xu)
01279 *     **************************************
01280       IMPLICIT NONE
01281       INCLUDE 'GLK.h'
01282       INTEGER           id,nchx
01283       DOUBLE PRECISION  xl,xu
01284       CHARACTER*80 title
01285       INTEGER           lact,ist,ist2
01286 *----------------------------------------------------------------------
01287       CALL GLK_hadres(id,lact)
01288       IF(lact .EQ. 0) THEN
01289          CALL GLK_Stop1('+++STOP in GLK_hinbo1: wrong id=',id)
01290       ENDIF
01291       ist=m_index(lact,2)
01292       ist2   = ist+7
01293       nchx   = m_b(ist2 +1)
01294       xl     = m_b(ist2 +2)
01295       xu     = m_b(ist2 +3)
01296       title  = m_titlc(lact)
01297       END
01298 
01299       SUBROUTINE GLK_UnPak(id,a,chd1,idum)
01300 *     *********************************
01301 * getting out histogram content (and error)
01302 * chd1= 'ERRO' is nonstandard option (unpack errors)
01303 *     ***********************************
01304       IMPLICIT NONE
01305       INCLUDE 'GLK.h'
01306       INTEGER           id,idum
01307       DOUBLE PRECISION  a(*)
01308       CHARACTER*(*) chd1
01309 *
01310       INTEGER            lact,ist,ist2,iflag2,ityphi,local,nch,nchy,ib
01311 *------------------------------------------------------------------------
01312       CALL GLK_hadres(id,lact)
01313       IF(lact .EQ. 0) goto 900
01314       ist   = m_index(lact,2)
01315       ist2  = ist+7
01316       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01317       ityphi   = mod(iflag2,10)
01318       IF(ityphi .EQ. 1) THEN
01319          nch   = m_b(ist2 +1)
01320          local = ist +m_buf1
01321       ELSEIF(ityphi .EQ. 2) THEN
01322          nchy  = m_b(ist2+5)
01323          nch   = nch*nchy
01324          local = ist+ m_buf2
01325       ELSE
01326          CALL GLK_Stop1('+++GLK_UnPak: check type of histo id=',id)
01327       ENDIF
01328       do 10 ib=1,nch
01329       IF(chd1 .NE. 'ERRO') THEN
01330 * normal bin
01331         a(ib) = m_b(local+ib)
01332       ELSE
01333 * error content
01334         IF(ityphi .EQ. 2) goto 901
01335         a(ib) = dsqrt( dabs(m_b(local+nch+ib) ))
01336       ENDIF
01337    10 continue
01338       RETURN
01339  900  CALL GLK_Retu1('+++GLK_UnPak: nonexisting id=',id)
01340       RETURN
01341  901  CALL GLK_Retu1('+++GLK_UnPak: no errors, two-dim, id=',id)
01342       END
01343 
01344       SUBROUTINE GLK_Pak(id,a)
01345 *     ************************
01346 * Loading in histogram content
01347 *     ***********************************
01348       IMPLICIT NONE
01349       INCLUDE 'GLK.h'
01350       INTEGER           id
01351       DOUBLE PRECISION  a(*)
01352 *
01353       INTEGER           lact,ist,ist2,iflag2,ityphi,nch,local,ib,nchy
01354 *----------------------------------------------------
01355       CALL GLK_hadres(id,lact)
01356       IF(lact .EQ. 0) goto 900
01357       ist  = m_index(lact,2)
01358       ist2 = ist+7
01359 * 2-dimens histo alowed
01360       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01361       ityphi   = mod(iflag2,10)
01362       IF(ityphi .EQ. 1) THEN
01363          nch   = m_b(ist2 +1)
01364          local = ist+m_buf1
01365       ELSEIF(ityphi .EQ. 2) THEN
01366          nchy  = m_b(ist2+5)
01367          nch   = nch*nchy
01368          local = ist+m_buf2
01369       ELSE
01370          CALL GLK_Stop1('+++GLK_Pak: wrong histo type, id=',id)
01371       ENDIF
01372       DO ib=1,nch
01373          m_b(local +ib) = a(ib)
01374       ENDDO
01375 * one nominal entry recorded
01376       m_index(lact,3)  = 1
01377       RETURN
01378   900 CONTINUE
01379       CALL GLK_Stop1('+++GLK_Pak: nonexisting id=',id)
01380       END
01381 
01382       SUBROUTINE GLK_Pake(id,a)
01383 *     **********************
01384 * Loading in error content
01385 *     ***********************************
01386       IMPLICIT NONE
01387       INCLUDE 'GLK.h'
01388       INTEGER           id
01389       DOUBLE PRECISION  a(*)
01390 *
01391       INTEGER           lact,ist,ist2,nch,iflag2,ityphi
01392       INTEGER           nb,ib
01393 *---------------------------------------------------------------------
01394       CALL GLK_hadres(id,lact)
01395       IF(lact .EQ. 0) goto 901
01396       ist  = m_index(lact,2)
01397       ist2 = ist+7
01398       nch=m_b(ist2+1)
01399 * 2-dimens histo NOT alowed
01400       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01401       ityphi   = mod(iflag2,10)
01402       IF(ityphi .NE. 1) GOTO 900
01403       DO ib=1,nch
01404          m_b(ist+m_buf1+nch+ib) = a(ib)**2
01405       ENDDO
01406       CALL GLK_idopt( id,'ERRO')
01407       RETURN
01408   900 CALL GLK_Stop1('+++GLK_Pake: only for 1-dim hist, id=',id)
01409       RETURN
01410   901 CALL GLK_Stop1('+++GLK_Pake: nonexisting id=',id)
01411       END
01412 
01413 
01414       SUBROUTINE GLK_Range1(id,ylr,yur)
01415 *     *****************************
01416 * provides y-scale for 1-dim plots
01417 *     ***********************************
01418       IMPLICIT NONE
01419       INCLUDE 'GLK.h'
01420       INTEGER           id
01421       DOUBLE PRECISION  ylr,yur
01422 *
01423       INTEGER           lact,ist,ist2,nch,ib,ioplog,iopsla,ioperb,iopsc1,iopsc2
01424       DOUBLE PRECISION  yl,yu
01425 *-------------------------------------------------------------
01426       CALL GLK_hadres(id,lact)
01427       IF(lact .EQ. 0) RETURN
01428       ist  = m_index(lact,2)
01429       ist2 = ist+7
01430       nch  = m_b(ist2 +1)
01431       yl   = m_b(ist+m_buf1+1)
01432       yu   = m_b(ist+m_buf1+1)
01433       DO ib=1,nch
01434          yl = min(yl,m_b(ist+m_buf1+ib))
01435          yu = max(yu,m_b(ist+m_buf1+ib))
01436       ENDDO
01437 * For default range some safety range is added
01438       yu = yu + 0.05*ABS(yu-yl)
01439 ***   yl = yl - 0.05*ABS(yu-yl) ! to be decided later on
01440 
01441 * If range was pre-defined then yl,yu are overwritten
01442       CALL GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
01443       IF(iopsc1 .EQ. 2) yl= m_b( ist +5)
01444       IF(iopsc2 .EQ. 2) yu= m_b( ist +6)
01445       ylr = yl
01446       yur = yu
01447       END
01448 
01449 
01450       SUBROUTINE GLK_hinbo2(id,nchx,xl,xu,nchy,yl,yu)
01451 *     *******************************************
01452       IMPLICIT NONE
01453       INCLUDE 'GLK.h'
01454       INTEGER           id,nchx,nchy
01455       DOUBLE PRECISION  xl,xu,yl,yu
01456       INTEGER           lact,ist,ist2
01457 *--------------------------------------------------
01458       CALL GLK_hadres(id,lact)
01459       IF(lact .EQ. 0) goto 900
01460       ist  = m_index(lact,2)
01461       ist2 = ist+7
01462       nchx = m_b(ist2 +1)
01463       xl   = m_b(ist2 +2)
01464       xu   = m_b(ist2 +3)
01465       nchy = m_b(ist2 +5)
01466       yl   = m_b(ist2 +6)
01467       yu   = m_b(ist2 +7)
01468       RETURN
01469   900 CALL GLK_Stop1(' +++GLK_hinbo2: nonexisting histo id= ',id)
01470       END
01471 
01472 
01473       SUBROUTINE GLK_Ymaxim(id,wmax)
01474 *     **************************
01475       IMPLICIT NONE
01476       INCLUDE 'GLK.h'
01477       INTEGER           id
01478       DOUBLE PRECISION  wmax
01479       INTEGER           lact,ist,jd,k
01480 *-------------------------------------------------------
01481       IF(id .NE. 0) THEN
01482          CALL GLK_hadres(id,lact)
01483          IF(lact .EQ. 0) RETURN
01484          ist= m_index(lact,2)
01485          m_b(ist+6) =wmax
01486          CALL GLK_idopt(id,'YMAX')
01487       ELSE
01488          DO k=1,m_idmax
01489             IF(m_index(k,1) .EQ. 0) GOTO 20
01490             ist=m_index(k,2)
01491             jd =m_index(k,1)
01492             m_b(ist+6) =wmax
01493             CALL GLK_idopt(jd,'YMAX')
01494          ENDDO
01495  20      CONTINUE
01496       ENDIF
01497       END
01498 
01499       SUBROUTINE GLK_Yminim(id,wmin)
01500 *     ******************************
01501       IMPLICIT NONE
01502       INCLUDE 'GLK.h'
01503       INTEGER           id
01504       DOUBLE PRECISION  wmin
01505       INTEGER           lact,ist,k,jd
01506 *---------------------------------------------
01507       IF(id .NE. 0) THEN
01508          CALL GLK_hadres(id,lact)
01509          IF(lact .EQ. 0) RETURN
01510          ist =m_index(lact,2)
01511          m_b(ist+5) =wmin
01512          CALL GLK_idopt(id,'YMIN')
01513       ELSE
01514          DO k=1,m_idmax
01515             IF(m_index(k,1) .EQ. 0) GOTO 20
01516             ist=m_index(k,2)
01517             jd =m_index(k,1)
01518             m_b(ist+5) =wmin
01519             CALL GLK_idopt(jd,'YMIN')
01520          ENDDO
01521  20      CONTINUE
01522       ENDIF
01523       END
01524 
01525       SUBROUTINE GLK_Reset(id,chd1)
01526 *     **************************
01527       IMPLICIT NONE
01528       INCLUDE 'GLK.h'
01529       INTEGER       id
01530       CHARACTER*(*) chd1
01531       INTEGER  lact,ist,ist2,iflag2,ityphi,ist3,nchx,nch,local,nchy,j
01532 *-------------------------------------------
01533       CALL GLK_hadres(id,lact)
01534       IF(lact .LE. 0) RETURN
01535       ist  =m_index(lact,2)
01536       ist2 = ist+7
01537 *
01538       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01539       ityphi   = mod(iflag2,10)
01540       IF(ityphi .EQ. 1) THEN
01541 * one-dim.
01542         ist3  = ist+11
01543         nchx  = m_b(ist2 +1)
01544         nch   = 2*nchx
01545         local = ist + m_buf1
01546       ELSEIF(ityphi .EQ. 2) THEN
01547 * two-dim.
01548         ist3  = ist+15
01549         nchx  = m_b(ist2 +1)
01550         nchy  = m_b(ist2 +5)
01551         nch   = nchx*nchy
01552         local = ist +m_buf2
01553       ELSE
01554          CALL GLK_Stop1('+++GLK_Reset: wrong type  id=',id)
01555       ENDIF
01556 * reset miscaelaneous entries and bins
01557       DO j=ist3+1,local +nch
01558          m_b(j)    = 0d0
01559       ENDDO
01560 * and no. of entries in m_index
01561       m_index(lact,3) = 0
01562       END
01563 
01564       SUBROUTINE GLK_Delet(id1)
01565 *     *********************
01566 * Now it should work (stj Nov. 91) but watch out!
01567 * should works for 2-dim histos, please check this!
01568 *     ***********************************
01569       IMPLICIT NONE
01570       INCLUDE 'GLK.h'
01571       INTEGER           id1
01572 *
01573       LOGICAL GLK_Exist
01574       INTEGER           id,lact,ist,ist2,nch,iflag2,ityphi,local,k,i,l,next,idec,nchx,nchy
01575 *--------------------------------------------
01576       ID=ID1
01577       IF(id .EQ. 0) GOTO 300
01578       IF( .NOT. GLK_Exist(id)) GOTO 900
01579       CALL GLK_hadres(id,lact)
01580       ist  = m_index(lact,2)
01581       ist2 = ist+7
01582 *----
01583 *[[[      WRITE(6,*) 'GLK_DELET-ing ID= ',ID
01584       idec    = nint(m_b(ist+2)-9d0-9d12)/10
01585       IF(idec .NE. id) WRITE(6,*)
01586      $     '++++GLK_DELET: ALARM! ID,IDEC= ',ID,IDEC
01587 *----
01588       nch  = m_b(ist2 +1)
01589       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
01590       ityphi   = MOD(iflag2,10)
01591       IF(ityphi .EQ. 1) THEN
01592 * one-dim.
01593         nchx  = m_b(ist2 +1)
01594         nch   = 2*nchx
01595 * lenght of local histo to be removed
01596         local = nch+m_buf1+1
01597       ELSEIF(ityphi .EQ. 2) THEN
01598 * two-dim.
01599         nchx  = m_b(ist2 +1)
01600         nchy  = m_b(ist2 +5)
01601         nch   = nchx*nchy
01602 * lenght of local histo to be removed
01603         local = nch+m_buf2+1
01604       ELSE
01605          CALL GLK_Stop1('+++GLK_Delet: wrong type id=',id)
01606       ENDIF
01607 * starting position of next histo in storage b
01608       next = ist+1 +local
01609 * move down all histos above this one
01610       DO 15 k =next,m_length
01611       m_b(k-local)=m_b(k)
01612    15 CONTINUE
01613 * define new end of storage
01614       m_length=m_length-local
01615 * clean free space at the end of storage b
01616       DO 20 k=m_length+1, m_length+local
01617    20 m_b(k)=0d0
01618 * shift adresses of all displaced histos
01619       DO 25 l=lact+1,m_idmax
01620       IF(m_index(l,1) .NE. 0) m_index(l,2)=m_index(l,2)-local
01621    25 CONTINUE
01622 * move entries in m_index down by one and remove id=lact entry
01623       DO 30 l=lact+1,m_idmax
01624       m_index(l-1,1)=m_index(l,1)
01625       m_index(l-1,2)=m_index(l,2)
01626       m_index(l-1,3)=m_index(l,3)
01627       m_titlc(l-1)=m_titlc(l)
01628    30 CONTINUE
01629 * last entry should be always empty
01630       m_index(m_idmax,1)=0
01631       m_index(m_idmax,2)=0
01632       m_index(m_idmax,3)=0
01633       do 50 k=1,80
01634    50 m_titlc(m_idmax)(k:k)=' '
01635       RETURN
01636 * -----------------------------------
01637 * Deleting all histos at once!!!
01638   300 m_length=0
01639       DO 400 i=1,m_idmax
01640       DO 340 k=1,3
01641   340 m_index(i,k)=0
01642       DO 350 k=1,80
01643   350 m_titlc(i)(k:k)=' '
01644  400  CONTINUE
01645       RETURN
01646 * -----------------------------------
01647  900  CONTINUE
01648       CALL GLK_Retu1(' +++GLK_DELET: nonexisting histo id= ',id)
01649       END
01650 
01651 
01652       SUBROUTINE GLK_Copch(ch1,ch2)
01653 *     *****************************
01654       IMPLICIT NONE
01655 * copies CHARACTER*80 ch1 into ch2 up to a first $ sign
01656       CHARACTER*80 ch1,ch2
01657       LOGICAL met
01658       INTEGER      i
01659 *----------------------------
01660       met = .FALSE.
01661       DO i=1,80
01662          IF( ch1(i:i) .EQ. '$' .or. met )   THEN
01663             ch2(i:i)=' '
01664             met=.TRUE.
01665          ELSE
01666             ch2(i:i)=ch1(i:i)
01667          ENDIF
01668       ENDDO
01669       END
01670 
01671       INTEGER FUNCTION GLK_jadre2(id)
01672 *------------------------------------------------
01673 * Good old version -- but it is very very slow!!!
01674 * In the case of 100 histograms or more.
01675 *------------------------------------------------
01676       IMPLICIT NONE 
01677       INCLUDE 'GLK.h'
01678       INTEGER           id,i
01679 *---------------------------------------
01680       GLK_jadre2=0
01681       DO 1 i=1,m_idmax
01682       IF(m_index(i,1) .EQ. id) goto 2
01683     1 CONTINUE
01684 * Nothing found.
01685       RETURN
01686 * Found: id=0 is also legitimate find!!!
01687     2 GLK_jadre2=i
01688       END
01689 
01690       SUBROUTINE GLK_hadres(id1,jadres)
01691 *     *********************************
01692 *--------------------------------------------------------------------
01693 * Educated guess based on past history is used to find quickly
01694 * location of the histogram in the matrix m_index.
01695 * This is based on observation that subsequent histogram calls
01696 * are linked into loops (so one can predict easily which histo will
01697 * be called next time).
01698 *--------------------------------------------------------------------
01699       IMPLICIT NONE
01700       INCLUDE 'GLK.h'
01701       INTEGER       id1,jadres
01702       INTEGER       ist,i,id
01703 *----------------------------------------------------------------------
01704       INTEGER iguess,jdlast,idlast
01705       DATA    iguess,jdlast,idlast /-2141593,-3141593,-3141593/
01706       SAVE    iguess,jdlast,idlast
01707 *----------------------------------------------------------------------
01708       id=id1
01709 * --- The case of ID=0 treated separately, it is used to find out
01710 * --- last entry in the m_index (it is marked with zero)
01711       IF(id .EQ. 0) THEN
01712          DO i=1,m_idmax
01713             IF(m_index(i,1) .EQ. 0) GOTO 44
01714          ENDDO
01715          CALL GLK_Stop1('+++GLK_hadres: Too short m_index=',m_index)
01716  44      CONTINUE
01717          jadres = i
01718          RETURN
01719       ENDIF
01720 
01721 * --- Omit sophistications if lack of initialization
01722       IF(jdlast .EQ. -3141593) GOTO 10
01723       IF(iguess .EQ. -2141593) GOTO 10
01724       IF(iguess .EQ. 0) GOTO 10
01725       IF(jdlast .EQ. 0) GOTO 10
01726 
01727 * --- Try first previous histo (for repeated calls)
01728       IF(jdlast .LT. 1 .OR. jdlast .GT. m_idmax) THEN
01729          WRITE(6,*) '+++++ hadres: jdlast=',jdlast
01730       ENDIF
01731       IF(m_index(jdlast,1) .EQ. id) THEN
01732          jadres = jdlast
01733 *##   WRITE(6,*)
01734 *##   $   'found, guess based on previous call to jadres ',jdlast
01735          GOTO 20
01736       ENDIF
01737 
01738 * --- Try current guess based on previous call
01739       IF(iguess .LT. 1 .OR. iguess .GT. m_idmax)  THEN
01740          WRITE(6,*)'+++++ hadres: iguess=',iguess
01741       ENDIF
01742       IF(m_index(iguess,1) .EQ. id) THEN
01743          jadres = iguess
01744 *##   WRITE(6,*)
01745 *##   $   'found, guess on previous calls recorded in m_b(ist+7)',jdlast
01746          GOTO 20
01747       ENDIF
01748 
01749 * ================================================
01750 *    Do it HARD WAY, Search all matrix m_index
01751 * ================================================
01752  10   CONTINUE
01753 *##   WRITE(6,*) 'trying HARD WAY'
01754       DO i=1,m_idmax
01755          jadres=i
01756          IF(m_index(i,1) .EQ. id) GOTO 20
01757       ENDDO
01758 * -------------------------------------
01759 *     Nothing found: jadres=0
01760 * -------------------------------------
01761       jadres=0
01762       RETURN
01763 * =====================================
01764 *     Found: Set new guess for next call
01765 * =====================================
01766  20   CONTINUE
01767 * --- and store result as a new guess in previous histo
01768 * --- but only if it existed!!!!
01769       DO i=1,m_idmax
01770          IF(m_index(i,1) .EQ. 0) GOTO 40
01771          IF(m_index(i,1) .EQ. idlast) THEN
01772             ist=m_index(i,2)
01773             IF(ist .GT. 0 .AND. ist .LT. m_LenmB) m_b(ist +7) = jadres
01774 *##   WRITE(6,*) 'STORED     id=',id
01775             GOTO 40
01776          ENDIF
01777       ENDDO
01778  40   CONTINUE
01779 *##   WRITE(6,*)  'found, hard way searching all of m_index)', jdlast
01780       iguess = m_b( m_index(jadres,2) +7)
01781       jdlast = jadres
01782       idlast = id
01783       END
01784 
01785 
01786 *--------------------------------------------------------------
01787 * ----------- storing histograms in the disk file -------------
01788 *--------------------------------------------------------------
01789 
01790       SUBROUTINE GLK_ReadFile(Dname)
01791 *     ******************************
01792       IMPLICIT NONE
01793       INCLUDE 'GLK.h'
01794       SAVE
01795       INTEGER ninph
01796       CHARACTER*60 Dname
01797 *-------------------------------------------------
01798       CALL GLK_Initialize
01799 * Read histograms
01800       WRITE(    *,*) 'GLK_ReadFile: Reading histos from ', Dname
01801       WRITE(m_out,*) 'GLK_ReadFile: Reading histos from ', Dname
01802       ninph = 21
01803       OPEN(ninph,file=Dname)           ! Open dump-file
01804       CALL GLK_hrfile(ninph,' ',' ')   ! Define dump-file unit in GKL
01805       CALL GLK_hrin(0,0,0)             ! Read histos from dump-file
01806       CALL GLK_hrend(' ')              ! Close dump-file
01807       END
01808 
01809       SUBROUTINE GLK_WriteFile(Dname)
01810 *     ******************************
01811       IMPLICIT NONE
01812       INCLUDE 'GLK.h'
01813       SAVE
01814       INTEGER nouth
01815       CHARACTER*60 Dname
01816 *-------------------------------------------------
01817       CALL GLK_Initialize
01818 * Write all histograms into disk file
01819       WRITE(    *,*) 'GLK_WriteFile: Writing histos into ', Dname
01820       WRITE(m_out,*) 'GLK_WriteFile: Writing histos into ', Dname
01821       nouth=22
01822       OPEN(nouth,file=Dname)           ! Open dump-file
01823       CALL GLK_hrfile(nouth,' ',' ')   ! Define dump-file in GLK
01824       CALL GLK_hrout(0,0,' ')          ! Dump all histos on disk
01825       CALL GLK_hrend(' ')              ! Close dump-file
01826       END
01827 
01828       SUBROUTINE GLK_hrfile(nhruni,chd1,chd2)
01829 *     ***************************************
01830       IMPLICIT NONE
01831       CHARACTER*(*) chd1,chd2
01832       INCLUDE 'GLK.h'
01833       SAVE
01834       INTEGER nhruni
01835 *-------------------------------------------------
01836       CALL GLK_Initialize
01837       m_huni=nhruni
01838       END
01839 
01840       SUBROUTINE GLK_hrout(idum1,idum2,chdum)
01841 *     ***************************************
01842       IMPLICIT NONE
01843       CHARACTER*8 chdum
01844 *
01845       INCLUDE 'GLK.h'
01846       SAVE
01847       INTEGER i,k,last,idum1,idum2
01848 *-----------------------------------------------------------------------
01849       CALL GLK_Initialize
01850       CALL GLK_hadres(0,last)
01851       IF(last.EQ.0) GOTO 900
01852       m_LenInd = last -1 
01853       WRITE(m_huni,'(6i10)')    m_version, m_LenInd, m_LenmB, m_Length
01854       WRITE(m_huni,'(6i10)')    ((m_index(i,k),k=1,3),i=1,m_LenInd)
01855       WRITE(m_huni,'(a80)')     (m_titlc(i),          i=1,m_LenInd)
01856       WRITE(m_huni,'(3d24.16)') (m_b(i), i=1,m_length)
01857       RETURN
01858  900  CONTINUE
01859       WRITE(m_out,*) '+++ GLK_hrout: no space in index'
01860       WRITE(    *,*) '+++ GLK_hrout: no space in index'
01861       END
01862 
01863 
01864       SUBROUTINE GLK_hrin(idum1,idum2,idum3)
01865 *     **************************************
01866 * New version which has a possibility to
01867 *            MERGE histograms
01868 * If given ID already exists then it is modified by adding 1000000 !!!!
01869 * Mergigng is done simply by appending new histograms at the
01870 * very end of the m_index and bin matrices.
01871 *     ***********************************
01872       IMPLICIT NONE
01873       INCLUDE 'GLK.h'
01874       INTEGER        idum1,idum2,idum3
01875       INTEGER        l,lact,lenold,istn,idn,k,lenind3,nvrs3,nouth
01876       INTEGER        i,lengt3,lenma3
01877 * Copy of the new m_index from the disk
01878       INTEGER        lndex(m_idmax,3)
01879       CHARACTER*80   titld(m_idmax)
01880       LOGICAL GLK_Exist
01881 *-----------------------------------------------------------
01882       CALL GLK_Initialize
01883       nouth=m_huni
01884 * Read basic params
01885       READ(nouth,'(6i10)')   nvrs3,lenind3,lenma3,lengt3
01886       IF(m_length+lengt3 .GE. m_LenmB) GOTO 900
01887 * Check version
01888       IF(m_version .NE. nvrs3) WRITE(m_out,*)
01889      $ '+++++WARNING (GLK_hrin): histos produced by older version',nvrs3
01890       IF(m_version .NE. nvrs3) WRITE(6,*)
01891      $ '+++++WARNING (GLK_hrin): histos produced by older version',nvrs3
01892       DO i=1,m_idmax
01893          DO k=1,3
01894             lndex(i,k)=0
01895          ENDDO
01896       ENDDO
01897 * We keep backward compatibility for reading disk files
01898       IF(nvrs3. LT. 130) lenind3 = m_idmax
01899 * Read new index and title from the disk
01900       READ(nouth,'(6i10)')  ((lndex(i,k),k=1,3),i=1,lenind3)
01901       READ(nouth,'(a80)')   (titld(i),          i=1,lenind3)
01902       lenold=m_length
01903 * Append AT ONCE content of new histos at the end of storage m_b
01904       m_length=m_length+lengt3
01905       READ(nouth,'(3d24.16)') (m_b(i),i=lenold+1,m_length)
01906 
01907 * Append ONE BY ONE m_index and m_titlc with new histos
01908       CALL GLK_hadres(0,lact)
01909       DO l=1,lenind3
01910          IF(lact .EQ. 0) GOTO 901
01911          idn= lndex(l,1)
01912          IF(idn .EQ. 0) GOTO 100
01913 * Identical id's are changed by adding big number = 1000000
01914  10      CONTINUE
01915          IF( GLK_Exist(idn) ) THEN
01916             idn = idn +1000000*(idn/iabs(idn))
01917             GOTO 10
01918          ENDIF
01919          m_index(lact,1)=idn
01920          m_index(lact,2)=lndex(l,2)+lenold
01921          m_index(lact,3)=lndex(l,3)
01922          m_titlc(lact)  =titld(l)
01923 * Still one small correction in the newly appended histo
01924          istn  = m_index(lact,2)
01925          m_b(istn +2)  = 9d12 +     idn*10 +9d0
01926          lact=lact+1
01927       ENDDO
01928   100 CONTINUE
01929       RETURN
01930 
01931  900  CONTINUE
01932       CALL GLK_Stop1('++++ GLK_hrin: to litle space, m_LenmB= ',m_LenmB)
01933  901  CONTINUE
01934       CALL GLK_Stop1('++++ GLK_hrin: to many histos, m_idmax= ',m_idmax)
01935       END
01936 
01937 
01938       SUBROUTINE GLK_hrin2(idum1,idum2,idum3)
01939 *     **********************************
01940 * New version which has a possibility to
01941 *            ADD histograms
01942 * If ID is not existing already then no action is taken
01943 *     ***********************************
01944       IMPLICIT NONE
01945       INCLUDE 'GLK.h'
01946       INTEGER         idum1,idum2,idum3
01947 * Copy of the histos from the disk
01948       DOUBLE PRECISION  bz(m_LenmB)
01949       INTEGER           indez(m_idmax,3)
01950       CHARACTER*80      titlz(m_idmax)
01951       LOGICAL           GLK_Exist
01952       INTEGER           nouth,ist3,nchx,ist,ist2,ist3z,nchxz,istz
01953       INTEGER           ist2z,lact,lenmaz,lengtz,nvrsz,lenindz,lz,id,i,k
01954 *-------------------------------------------------------------
01955       CALL GLK_Initialize
01956       nouth=m_huni
01957 * Read basic params
01958       READ(nouth,'(6i10)')   nvrsz,lenindz,lenmaz,lengtz
01959 * Check version
01960       IF(m_version .NE. nvrsz) WRITE(m_out,*)
01961      $ '++++WARNING (GLK_hrin2): histos produced by older version',nvrsz
01962       IF(m_version .NE. nvrsz) WRITE(6,*)
01963      $ '++++WARNING (GLK_hrin2): histos produced by older version',nvrsz
01964 
01965 * We keep backward compatibility for reading disk files
01966       IF(nvrsz. LT. 130) lenindz = m_idmax
01967       DO i=1,m_idmax
01968          DO k=1,3
01969             indez(i,k)=0
01970          ENDDO
01971       ENDDO
01972 * Read new m_index, title and bins from the disk
01973       READ(nouth,'(6i10)')    ((indez(i,k),k=1,3),i=1,lenindz)
01974       READ(nouth,'(a80)')     (titlz(i) ,         i=1,lenindz)
01975       READ(nouth,'(3d24.16)') (bz(i),i=1,lengtz)
01976 
01977 * Add new histos from disk to existing ones one by one
01978       DO lz=1,lenindz
01979          id= indez(lz,1)
01980          IF(id .EQ. 0) GOTO 200
01981          IF( .NOT. GLK_Exist(id)) THEN
01982             CALL GLK_Retu1('GLK_hrin2: unmached skipped histo ID=', id)
01983             GOTO 100
01984          ENDIF
01985 * parameters of existing histo
01986          CALL GLK_hadres(id,lact)
01987          ist  = m_index(lact,2)
01988          ist2 = ist+7
01989          ist3 = ist+11
01990          nchx = m_b(ist2 +1)
01991 * parameters of the histo from the disk
01992          istz   = indez(lz,2)
01993          ist2z  = istz+7
01994          ist3z  = istz+11
01995          nchxz  = bz(ist2z +1)
01996          IF(nchx .NE. nchxz) THEN
01997             CALL GLK_Retu1('GLK_hrin2: non-equal bins, skipped ID=', id)
01998             GOTO 100
01999          ENDIF
02000 * Add/Merge all additive entries of the two histos
02001 * No of entries in m_index
02002          m_index(lact,3) = m_index(lact,3)+indez(lact,3)
02003 * Overflows, underflows etc.
02004          DO i=1,12
02005             m_b(ist3+i)=m_b(ist3+i) +bz(ist3z+i)
02006          ENDDO
02007 * Except of this one non-additive entry
02008          m_b(ist3+13)=max(m_b(ist3+13),m_b(ist3z+13))
02009 * Regular bin content added now!
02010          DO i= 1, 2*nchx
02011             m_b(ist+m_buf1+i)=m_b(ist+m_buf1+i) +bz(istz+m_buf1+i)
02012          ENDDO
02013  100     CONTINUE
02014       ENDDO
02015  200  CONTINUE
02016       END
02017 
02018       SUBROUTINE GLK_hrend(chdum)
02019 *     ***********************
02020       IMPLICIT NONE
02021       INCLUDE 'GLK.h'
02022       CHARACTER*(*) chdum
02023 *---------------------------
02024       CLOSE(m_huni)
02025       END
02026 *======================================================================
02027 *                End of histograming procedures
02028 *======================================================================
02029 
02030 
02031 
02032 *======================================================================
02033 *               Ploting procedures using LaTeX
02034 *======================================================================
02035 
02036       SUBROUTINE GLK_PlInitialize(Lint,TeXfile)
02037 *----------------------------------------------------------------------
02038 * Lint =0 
02039 *     Normal mode, full LaTeX header
02040 * Lint =1
02041 *     For TeX file is used in \input, no  LaTeX header
02042 * Lint =2
02043 *     LaTeX header for one-page plot used as input for postscript
02044 *
02045 * Negative Lint only for debug, big frame around plot is added.
02046 *----------------------------------------------------------------------
02047       IMPLICIT NONE
02048       INCLUDE 'GLK.h'
02049       SAVE
02050       INTEGER Lint,noufig
02051       CHARACTER*60  TeXfile
02052 *---------------------------------
02053 * Initialize GLK_Plot
02054       CALL GLK_PlInt(Lint)             ! Define header style
02055       noufig=11                        ! 
02056       OPEN(noufig,file=TeXfile)        ! Open LaTeX file
02057       CALL GLK_PlCap(noufig)           ! Initialize GLK_Plot
02058       END
02059 
02060       SUBROUTINE GLK_PlEnd
02061 *     ********************
02062       IMPLICIT NONE
02063       INCLUDE 'GLK.h'
02064       SAVE
02065 *---------------------------------------------------
02066 * Note that TeX file is used in \input then you may not want
02067 * to have header and \end{document}
02068       IF( ABS(m_lint) .NE. 1) THEN
02069          WRITE(m_ltx,'(2A)') m_BS,'end{document}'
02070       ENDIF
02071       CLOSE(m_ltx)
02072       END
02073 
02074       SUBROUTINE GLK_PlInt(Lint)
02075 *     **************************
02076       IMPLICIT NONE
02077       INCLUDE 'GLK.h'
02078       SAVE
02079       INTEGER Lint
02080 *---------------------------------
02081       m_lint = Lint
02082       END
02083 
02084       SUBROUTINE GLK_PlCap(LtxUnit)
02085 *     ****************************
02086       IMPLICIT NONE
02087       INCLUDE 'GLK.h'
02088       INTEGER   LtxUnit,i,k
02089 *---------
02090       CALL GLK_Initialize
02091       m_KeyTit= 0
02092       DO i=1,m_titlen
02093          DO k=1,80
02094             m_titch(i)(k:k)=' '
02095          ENDDO
02096       ENDDO
02097 *---------
02098       m_tline = 1
02099       m_ltx=IABS(LtxUnit)
02100 
02101       IF( ABS(m_lint) .EQ. 0) THEN
02102 * Normal mode, no colors!!!
02103          WRITE(m_ltx,'(A,A)') m_BS,'documentclass[12pt]{article}'
02104 *!!         WRITE(m_ltx,'(A,A)') m_BS,'usepackage{html}'
02105          WRITE(m_ltx,'(A,A)') m_BS,'textwidth  = 16cm'
02106          WRITE(m_ltx,'(A,A)') m_BS,'textheight = 18cm'
02107          WRITE(m_ltx,'(A,A)') m_BS,'begin{document}'
02108          WRITE(m_ltx,'(A)') '  '
02109       ELSEIF( ABS(m_lint) .EQ. 1) THEN
02110 * For TeX file is used in \input
02111          WRITE(m_ltx,'(A)') '  '
02112       ELSEIF( ABS(m_lint) .EQ. 2) THEN
02113 * For one-page plot being input for postrscript
02114 *!!         WRITE(m_ltx,'(A,A)') m_BS,'documentclass[12pt,html]{article}'
02115 *!!         WRITE(m_ltx,'(A,A)') m_BS,'documentclass[12pt,dvips]{seminar}' !<-for colors!!!
02116          WRITE(m_ltx,'(A,A)') m_BS,'documentclass[12pt,dvips]{article}'
02117          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{amsmath}'
02118          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{amssymb}'
02119 *!!         WRITE(m_ltx,'(A,A)') m_BS,'usepackage{html}'
02120          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{epsfig}'
02121          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{epic}'
02122          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{eepic}'
02123          WRITE(m_ltx,'(A,A)') m_BS,'usepackage{color}' !<-for colors!!!
02124 *!!         WRITE(m_ltx,'(A,A)') m_BS,'hoffset    = -1in'
02125 *!!         WRITE(m_ltx,'(A,A)') m_BS,'voffset    = -1in'
02126 *!!         WRITE(m_ltx,'(A,A)') m_BS,'textwidth  = 16cm'
02127 *!!         WRITE(m_ltx,'(A,A)') m_BS,'textheight = 16cm'
02128 *!!         WRITE(m_ltx,'(A,A)') m_BS,'oddsidemargin = 0cm'
02129 *!!         WRITE(m_ltx,'(A,A)') m_BS,'topmargin     = 0cm'
02130 *!!         WRITE(m_ltx,'(A,A)') m_BS,'headheight    = 0cm'
02131 *!!         WRITE(m_ltx,'(A,A)') m_BS,'headsep       = 0cm'
02132          WRITE(m_ltx,'(A,A)') m_BS,'begin{document}'
02133          WRITE(m_ltx,'(A,A)') m_BS,'pagestyle{empty}'
02134          WRITE(m_ltx,'(A)') '  '
02135       ELSE
02136          CALL GLK_Stop1('+++STOP in GLK_PlInt, wrong m_lint=',m_lint)
02137       ENDIF
02138       END
02139 
02140 
02141       SUBROUTINE GLK_Plot(id,ch1,ch2,kdum)
02142 *     ************************************
02143       IMPLICIT NONE
02144       INCLUDE 'GLK.h'
02145       CHARACTER CH1,CH2,CHR
02146       CHARACTER*80 TITLE
02147       INTEGER          id,kdum
02148       DOUBLE PRECISION YY(m_MaxNb),YER(m_MaxNb)
02149       LOGICAL GLK_EXIST
02150       INTEGER          idum,kax,kay,ioplog,iopsla,ioperb,iopsc1,iopsc2
02151       INTEGER          ker,nchx
02152       DOUBLE PRECISION XL,XU,DXL,DXU,yl,yu
02153 *--------------------------------------------
02154       DATA CHR /' '/
02155 * RETURN if histo non-existing
02156       IF(.NOT.GLK_EXIST(ID)) GOTO 900
02157 * ...unpack histogram
02158       CALL GLK_UNPAK(ID,YY ,'    ',IDUM)
02159       CALL GLK_UNPAK(ID,YER,'ERRO',IDUM)
02160       CALL GLK_HINBO1(ID,TITLE,NCHX,DXL,DXU)
02161       XL = DXL
02162       XU = DXU
02163       CALL GLK_RANGE1(ID,YL,YU)
02164       kax=1200
02165       kay=1200
02166       IF(CH1 .EQ. 'S') THEN
02167 * ...superimpose plot
02168         BACKSPACE(m_ltx)
02169         BACKSPACE(m_ltx)
02170       ELSE
02171 * ...new frame only
02172         CHR=CH1
02173         CALL GLK_Plfram1(ID,kax,kay)
02174       ENDIF
02175       WRITE(m_ltx,'(A)')     '%========== next plot (line) =========='
02176       WRITE(m_ltx,'(A,I10)') '%==== HISTOGRAM ID=',ID
02177       WRITE(m_ltx,'(A,A70 )')'% ',TITLE
02178 *...cont. line for functions
02179       CALL GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
02180       ker = ioperb-1
02181       IF (iopsla .EQ. 2)  CHR='C'
02182 *...suppress GLK_PLOT assignments
02183       IF (CH2 .EQ. 'B')   CHR=' '
02184       IF (CH2 .EQ. '*')   CHR='*'
02185       IF (CH2 .EQ. 'C')   CHR='C'
02186 *...various types of lines
02187       IF     (CHR .EQ. ' ') THEN
02188 *...contour line used for histogram
02189           CALL GLK_PlHist(kax,kay,NCHX,YL,YU,YY,KER,YER)
02190       ELSE IF(CHR .EQ. '*') THEN
02191 *...marks in the midle of the bin
02192           CALL GLK_PlHis2(kax,kay,NCHX,YL,YU,YY,KER,YER)
02193       ELSE IF(CHR .EQ. 'C') THEN
02194 *...slanted (dotted) line in plotting non-MC functions
02195           CALL GLK_PlCirc(kax,kay,NCHX,YL,YU,YY)
02196       ENDIF
02197 *------------------------------!
02198 * Ending
02199 *------------------------------!
02200       WRITE(m_ltx,'(2A)') m_BS,'end{picture} % close entire picture '
02201       WRITE(m_ltx,'(2A)') m_BS,'end{figure}'
02202 
02203       RETURN
02204   900 CALL GLK_Retu1('+++GLK_PLOT: Nonexistig histo, skipped, id=' ,ID)
02205       END
02206 
02207       SUBROUTINE GLK_Plfram1(ID,kax,kay)
02208 *     **********************************
02209       IMPLICIT NONE
02210       INCLUDE 'GLK.h'
02211       INTEGER           ID,kax,kay
02212       CHARACTER*80 title
02213       DOUBLE PRECISION   TIPSY(20),TIPSX(20)
02214       DOUBLE PRECISION   XL,DXL,XU,DXU
02215       INTEGER            ntipy,ntipx,nchx,icont
02216       DOUBLE PRECISION   yu,yl
02217       DATA ICONT/0/
02218 *----------------
02219       ICONT=ICONT+1
02220       CALL GLK_HINBO1(ID,TITLE,NCHX,DXL,DXU)
02221       XL = DXL
02222       XU = DXU
02223       CALL GLK_RANGE1(ID,YL,YU)
02224 
02225       IF(ICONT .GT. 1) WRITE(m_ltx,'(2A)') m_BS,'newpage'
02226 *------------------------------!
02227 *           Header
02228 *------------------------------!
02229       WRITE(m_ltx,'(A)') ' '
02230       WRITE(m_ltx,'(A)') ' '
02231       WRITE(m_ltx,'(A)') 
02232 '%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     $%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%'
02233       WRITE(m_ltx,'(A)') 
02234 '%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     $%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%'
02235       WRITE(m_ltx,'(2A)') m_BS,'begin{figure}[!ht]'
02236       WRITE(m_ltx,'(2A)') m_BS,'centering'
02237 *------------------------------!
02238 * General Caption
02239 *------------------------------!
02240       WRITE(m_ltx,'(4A)') m_BS,'caption{',m_BS,'small'
02241       IF(M_KEYTIT.EQ.0) THEN
02242         WRITE(m_ltx,'(A)')     TITLE
02243       ELSE
02244         WRITE(m_ltx,'(A)')     m_titch(1)
02245       ENDIF
02246       WRITE(m_ltx,'(A)') '}'
02247 *------------------------------!
02248 * Frames and labels
02249 *------------------------------!
02250       WRITE(m_ltx,'(A)') '% =========== big frame, title etc. ======='
02251       WRITE(m_ltx,'(4A)') m_BS,'setlength{',m_BS,'unitlength}{0.1mm}'
02252       WRITE(m_ltx,'(2A)') m_BS,'begin{picture}(1600,1500)'
02253       WRITE(m_ltx,'(4A)')
02254      $     m_BS,'put(0,0){',m_BS,'framebox(1600,1500){ }}'
02255       WRITE(m_ltx,'(A)') '% =========== small frame, labeled axis ==='
02256       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02257      $    m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
02258       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02259      $    m_BS,'put(0,0){',m_BS,'framebox( ',kax,',',kay,'){ }}'
02260       WRITE(m_ltx,'(A)') '% =========== x and y axis ================'
02261       CALL GLK_SAxisX(kax,XL,XU,NTIPX,TIPSX)
02262       CALL GLK_SAxisY(kay,YL,YU,NTIPY,TIPSY)
02263       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}'
02264      $                ,'% end of plotting labeled axis'
02265       END
02266 
02267       SUBROUTINE GLK_SAxisX(kay,YL,YU,NLT,TIPSY)
02268 *     ***************************************
02269 * plotting x-axis with long and short tips
02270       IMPLICIT NONE
02271       INCLUDE 'GLK.h'
02272       INTEGER           kay,NLT
02273       DOUBLE PRECISION  YL,YU,TIPSY(20)
02274 *
02275       INTEGER           LY,JY,n,nts,k,lex
02276       DOUBLE PRECISION  DY,pds,scmx,p0s,ddys,yy0l,ddyl,pdl,p0l,yy0s
02277 *---------------------------------------------------
02278       DY= ABS(YU-YL)
02279       LY = NINT( LOG10(DY) -0.4999999d0 )
02280       JY = NINT(DY/10d0**LY)
02281       DDYL = DY*10d0**(-LY)
02282       IF( JY .EQ. 1)             DDYL = 10d0**LY*0.25d0
02283       IF( JY .GE. 2.AND.JY .LE. 3) DDYL = 10d0**LY*0.5d0
02284       IF( JY .GE. 4.AND.JY .LE. 6) DDYL = 10d0**LY*1.0d0
02285       IF( JY .GE. 7)             DDYL = 10d0**LY*2.0d0
02286       WRITE(m_ltx,'(A)') '% .......GLK_SAxisX........ '
02287       WRITE(m_ltx,'(A,I4)') '%  JY= ',JY
02288 *-------
02289       NLT = INT(DY/DDYL)
02290       NLT = MAX0(MIN0(NLT,20),1)+1
02291       YY0L = NINT(YL/DDYL+0.5d0)*DDYL
02292       DDYS = DDYL/10d0
02293       YY0S = NINT(YL/DDYS+0.4999999d0)*DDYS
02294       P0L = kay*(YY0L-YL)/(YU-YL)
02295       PDL = kay*DDYL/(YU-YL)
02296       P0S = kay*(YY0S-YL)/(YU-YL)
02297       PDS = kay*DDYS/(YU-YL)
02298       NLT = INT(ABS(YU-YY0L)/DDYL+0.0000001d0)+1
02299       NTS = INT(ABS(YU-YY0S)/DDYS+0.0000001d0)+1
02300       DO 41 N=1,NLT
02301       TIPSY(N) =YY0L+ DDYL*(N-1)
02302   41  CONTINUE
02303       WRITE(m_ltx,1000)
02304      $ m_BS,'multiput('  ,P0L,  ',0)('  ,PDL,  ',0){'  ,NLT,  '}{',
02305      $ m_BS,'line(0,1){25}}',
02306      $ m_BS,'multiput('  ,P0S,  ',0)('  ,PDS,  ',0){'  ,NTS,  '}{',
02307      $ m_BS,'line(0,1){10}}'
02308       WRITE(m_ltx,1001)
02309      $ m_BS,'multiput('  ,P0L,  ','  ,kay,  ')('  ,PDL,  ',0){'  ,NLT,
02310      $ '}{'  ,m_BS,  'line(0,-1){25}}',
02311      $ m_BS,'multiput('  ,P0S,  ','  ,kay,  ')('  ,PDS,  ',0){'  ,NTS,
02312      $ '}{'  ,m_BS,  'line(0,-1){10}}'
02313  1000 FORMAT(2A,F8.2,A,F8.2,A,I4,3A)
02314  1001 FORMAT(2A,F8.2,A,I4,A,F8.2,A,I4,3A)
02315 * ...labeling of axis
02316       SCMX = DMAX1(DABS(YL),DABS(YU))
02317       LEX  = NINT( LOG10(SCMX) -0.50001)
02318       DO 45 N=1,NLT
02319       K = NINT(kay*(TIPSY(N)-YL)/(YU-YL))
02320       IF(LEX .LT. 2.AND.LEX .GT. -1) THEN
02321 * ...without exponent
02322       WRITE(m_ltx,'(2A,I4,5A,F8.3,A)')
02323      $ m_BS,'put(',K,',-25){',m_BS,'makebox(0,0)[t]{',m_BS,'large $ ',
02324      $ TIPSY(N), ' $}}'
02325       ELSE
02326 * ...with exponent
02327       WRITE(m_ltx,'(2A,I4,5A,F8.3,2A,I4,A)')
02328      $ m_BS,'put(' ,K, ',-25){',m_BS,'makebox(0,0)[t]{',m_BS,'large $ ',
02329      $ TIPSY(N)/(10d0**LEX),m_BS,'cdot 10^{',LEX,'} $}}'
02330       ENDIF
02331   45  CONTINUE
02332       END
02333 
02334       SUBROUTINE GLK_SAxisY(kay,yl,yu,nlt,tipsy)
02335 *     ******************************************
02336 * plotting y-axis with long and short tips
02337       IMPLICIT NONE
02338       INCLUDE 'GLK.h'
02339       INTEGER          kay,nlt
02340       DOUBLE PRECISION yl,yu,tipsy(20)
02341 *
02342       DOUBLE PRECISION dy,ddyl,z0l,scmx,yy0s,ddys,yy0l,p0l,pds,p0s,pdl
02343       INTEGER          ly,jy,n,nts,k,lex
02344 *---------------------------------------------------
02345       DY= ABS(YU-YL)
02346       LY = NINT( LOG10(DY) -0.49999999d0 )
02347       JY = NINT(DY/10d0**LY)
02348       DDYL = DY*10d0**(-LY)
02349       IF( JY .EQ. 1)             DDYL = 10d0**LY*0.25d0
02350       IF( JY .GE. 2.AND.JY .LE. 3) DDYL = 10d0**LY*0.5d0
02351       IF( JY .GE. 4.AND.JY .LE. 6) DDYL = 10d0**LY*1.0d0
02352       IF( JY .GE. 7)             DDYL = 10d0**LY*2.0d0
02353       WRITE(m_ltx,'(A)') '% .......GLK_SAxisY........ '
02354       WRITE(m_ltx,'(A,I4)') '%  JY= ',JY
02355 *-------
02356       NLT = INT(DY/DDYL)
02357       NLT = MAX0(MIN0(NLT,20),1)+1
02358       YY0L = NINT(YL/DDYL+0.4999999d0)*DDYL
02359       DDYS = DDYL/10d0
02360       YY0S = NINT(YL/DDYS+0.5d0)*DDYS
02361       P0L = kay*(YY0L-YL)/(YU-YL)
02362       PDL = kay*DDYL/(YU-YL)
02363       P0S = kay*(YY0S-YL)/(YU-YL)
02364       PDS = kay*DDYS/(YU-YL)
02365       NLT= INT(ABS(YU-YY0L)/DDYL+0.0000001d0) +1
02366       NTS= INT(ABS(YU-YY0S)/DDYS+0.0000001d0) +1
02367       DO 41 N=1,NLT
02368       TIPSY(N) =YY0L+ DDYL*(N-1)
02369   41  CONTINUE
02370 * plotting tics on vertical axis
02371       WRITE(m_ltx,1000)
02372      $ m_BS,'multiput(0,'  ,P0L,  ')(0,'  ,PDL  ,'){'  ,NLT,  '}{',
02373      $ m_BS,'line(1,0){25}}',
02374      $ m_BS,'multiput(0,'  ,P0S,  ')(0,'  ,PDS,  '){'  ,NTS,  '}{',
02375      $ m_BS,'line(1,0){10}}'
02376       WRITE(m_ltx,1001)
02377      $ m_BS,'multiput('  ,kay,  ','  ,P0L,  ')(0,'  ,PDL,  '){'  ,NLT,
02378      $ '}{',m_BS,'line(-1,0){25}}',
02379      $ m_BS,'multiput('  ,kay,  ','  ,P0S,  ')(0,'  ,PDS,  '){'  ,NTS,
02380      $ '}{',m_BS,'line(-1,0){10}}'
02381  1000 FORMAT(2A,F8.2,A,F8.2,A,I4,3A)
02382  1001 FORMAT(2A,I4,A,F8.2,A,F8.2,A,I4,3A)
02383 * ...Zero line if necessary
02384       Z0L = kay*(-YL)/(YU-YL)
02385       IF(Z0L .GT. 0D0.AND.Z0L .LT. FLOAT(kay))
02386      $      WRITE(m_ltx,'(2A,F8.2,3A,I4,A)')
02387      $       m_BS,'put(0,'  ,Z0L,  '){',m_BS,'line(1,0){'  ,kay,  '}}'
02388 * ...labeling of axis
02389       SCMX = DMAX1(DABS(YL),DABS(YU))
02390       LEX  = NINT( LOG10(SCMX) -0.50001d0)
02391       DO 45 N=1,NLT
02392       K = NINT(kay*(TIPSY(N)-YL)/(YU-YL))
02393       IF(LEX .LT. 2.AND.LEX .GT. -1) THEN
02394 * ...without exponent
02395       WRITE(m_ltx,'(2A,I4,5A,F8.3,A)')
02396      $  m_BS,'put(-25,'  ,K,  '){',m_BS,'makebox(0,0)[r]{',
02397      $  m_BS,'large $ '  ,TIPSY(N),  ' $}}'
02398       ELSE
02399 * ...with exponent
02400       WRITE(m_ltx,'(2A,I4,5A,F8.3,2A,I4,A)')
02401      $ m_BS,'put(-25,'  ,K,  '){',m_BS,'makebox(0,0)[r]{',
02402      $ m_BS,'large $ '
02403      $ ,TIPSY(N)/(10d0**LEX),  m_BS,'cdot 10^{'  ,LEX,  '} $}}'
02404       ENDIF
02405   45  CONTINUE
02406       END
02407 
02408       SUBROUTINE GLK_PlHist(kax,kay,nchx,yl,yu,yy,ker,yer)
02409 *     ************************************************
02410 * plotting contour line for histogram
02411 *     ***********************
02412       IMPLICIT NONE
02413       INCLUDE 'GLK.h'
02414       INTEGER           kax,kay,nchx,ker
02415       DOUBLE PRECISION  yl,yu,yy(*),yer(*)
02416       CHARACTER*80 FMT1
02417 *
02418       INTEGER           IX0,ix2,idx,ie,ierr,idy,ib,iy0,iy1,ix1
02419 *---------------------------------------------------
02420       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02421      $  m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
02422       WRITE(m_ltx,'(A)') '% ========== plotting primitives =========='
02423 *...various types of line
02424       IF(m_tline .EQ. 1) THEN
02425          WRITE(m_ltx,'(2A)') m_BS,'thicklines '
02426       ELSE
02427          WRITE(m_ltx,'(2A)') m_BS,'thinlines '
02428       ENDIF
02429 *...short macros for vertical/horizontal straight lines
02430       WRITE(m_ltx,'(8A)')
02431      $ m_BS,'newcommand{',m_BS,'x}[3]{',m_BS,'put(#1,#2){',
02432      $ m_BS,'line(1,0){#3}}}'
02433       WRITE(m_ltx,'(8A)')
02434      $ m_BS,'newcommand{',m_BS,'y}[3]{',m_BS,'put(#1,#2){',
02435      $ m_BS,'line(0,1){#3}}}'
02436       WRITE(m_ltx,'(8A)')
02437      $ m_BS,'newcommand{',m_BS,'z}[3]{',m_BS,'put(#1,#2){',
02438      $ m_BS,'line(0,-1){#3}}}'
02439 *   error bars
02440       WRITE(m_ltx,'(8A)')
02441      $   m_BS,'newcommand{',m_BS,'e}[3]{',
02442      $   m_BS,'put(#1,#2){',m_BS,'line(0,1){#3}}}'
02443       IX0=0
02444       IY0=0
02445       DO 100 IB=1,NCHX
02446       IX1 = NINT(kax*(IB-0.00001)/NCHX)   !ib=7
02447       IY1 = NINT(kay*(YY(IB)-YL)/(YU-YL)) !iy1=775,while ix0=168,iy0=770
02448       IDY = IY1-IY0
02449       IDX = IX1-IX0
02450       FMT1 = '(2(2A,I4,A,I4,A,I4,A))'
02451       IF( IDY .GE. 0) THEN
02452          IF(IY1 .GE. 0.AND.IY1 .LE. kay)
02453      $   WRITE(m_ltx,FMT1) m_BS,'y{',IX0,'}{',IY0,'}{',IDY,'}',
02454      $                     m_BS,'x{',IX0,'}{',IY1,'}{',IDX,'}'
02455       ELSE
02456          IF(IY1 .GE. 0.AND.IY1 .LE. kay)
02457      $   WRITE(m_ltx,FMT1) m_BS,'z{',IX0,'}{',IY0,'}{',-IDY,'}',
02458      $                     m_BS,'x{',IX0,'}{',IY1,'}{',IDX,'}'
02459       ENDIF
02460       IX0=IX1
02461       IY0=IY1
02462       IF(KER .EQ. 1) THEN
02463         IX2  = NINT(kax*(IB-0.5000d0)/NCHX)
02464         IERR = NINT(kay*((YY(IB)-YER(IB))-YL)/(YU-YL))
02465         IE = NINT(kay*YER(IB)/(YU-YL))
02466         IF(IY1 .GE. 0.AND.IY1 .LE. kay.and.abs(ierr) .LE. 9999
02467      $     .and.2*ie .LE. 9999) WRITE(m_ltx,8000) m_BS,IX2,IERR,IE*2
02468       ENDIF
02469  100  CONTINUE
02470 8000  FORMAT(4(A1,2He{,I4,2H}{,I5,2H}{,I4,1H}:1X ))
02471       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}',
02472      $       ' % end of plotting histogram'
02473 * change line-style
02474       m_tline= m_tline+1
02475       IF(m_tline .GT. 2) m_tline=1
02476       END
02477 
02478       SUBROUTINE GLK_PlHis2(kax,kay,nchx,yl,yu,yy,ker,yer)
02479 *     ************************************************
02480 * marks in the midle of the bin
02481 *     **********************************
02482       IMPLICIT NONE
02483       INCLUDE 'GLK.h'
02484       DOUBLE PRECISION  yl,yu,yy(*),yer(*)
02485       INTEGER           kax,kay,nchx,ker
02486 *
02487       INTEGER           iy1,ierr,ie,ix1,irad1,irad2,ib
02488 *---------------------------------------------------
02489 
02490       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02491      $ m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
02492       WRITE(m_ltx,'(A)') '% ========== plotting primitives =========='
02493 *...various types of mark
02494       IRAD1= 6
02495       IRAD2=10
02496       IF(m_tline .EQ. 1) THEN
02497 *   small filled circle
02498        WRITE(m_ltx,'(8A,I3,A)')
02499      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02500      $   m_BS,'put(#1,#2){',m_BS,'circle*{',IRAD1,'}}}'
02501       ELSEIF(m_tline .EQ. 2) THEN
02502 *   small open circle
02503        WRITE(m_ltx,'(8A,I3,A)')
02504      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02505      $   m_BS,'put(#1,#2){',m_BS,'circle{',IRAD1,'}}}'
02506       ELSEIF(m_tline .EQ. 3) THEN
02507 *   big filled circle
02508        WRITE(m_ltx,'(8A,I3,A)')
02509      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02510      $   m_BS,'put(#1,#2){',m_BS,'circle*{',IRAD2,'}}}'
02511       ELSEIF(m_tline .EQ. 4) THEN
02512 *   big open circle
02513        WRITE(m_ltx,'(8A,I3,A)')
02514      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02515      $   m_BS,'put(#1,#2){',m_BS,'circle{',IRAD2,'}}}'
02516 * Other symbols
02517       ELSEIF(m_tline .EQ. 5) THEN
02518        WRITE(m_ltx,'(10A)')
02519      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02520      $   m_BS,'put(#1,#2){',m_BS,'makebox(0,0){$',m_BS,'diamond$}}}'
02521       ELSE
02522        WRITE(m_ltx,'(10A)')
02523      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02524      $   m_BS,'put(#1,#2){',m_BS,'makebox(0,0){$',m_BS,'star$}}}'
02525       ENDIF
02526 *   error bars
02527       WRITE(m_ltx,'(8A)')
02528      $   m_BS,'newcommand{',m_BS,'E}[3]{',
02529      $   m_BS,'put(#1,#2){',m_BS,'line(0,1){#3}}}'
02530       DO 100 IB=1,NCHX
02531       IX1 = NINT(kax*(IB-0.5000d0)/NCHX)
02532       IY1 = NINT(kay*(YY(IB)-YL)/(YU-YL))
02533       IF(IY1 .GE. 0.AND.IY1 .LE. kay) WRITE(m_ltx,7000) m_BS,IX1,IY1
02534       IF(KER .EQ. 1) THEN
02535         IERR = NINT(kay*((YY(IB)-YER(IB))-YL)/(YU-YL))
02536         IE   = NINT(kay*YER(IB)/(YU-YL))
02537         IF(IY1 .GE. 0.AND.IY1 .LE. kay.and.abs(ierr) .LE. 9999
02538      $       .and.2*ie .LE. 9999) WRITE(m_ltx,8000) m_BS,IX1,IERR,IE*2
02539       ENDIF
02540  100  CONTINUE
02541 7000  FORMAT(4(A1,2HR{,I4,2H}{,I4,1H}:1X ))
02542 8000  FORMAT(4(A1,2HE{,I4,2H}{,I5,2H}{,I4,1H}:1X ))
02543       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}',
02544      $    ' % end of plotting histogram'
02545 * change line-style
02546       m_tline= m_tline+1
02547       IF(m_tline .GT. 6) m_tline=1
02548       END
02549 
02550       SUBROUTINE GLK_PlCirc(kax,kay,nchx,yl,yu,yy)
02551 *     ****************************************
02552 * plots equidistant points, four-point interpolation,
02553       IMPLICIT NONE
02554       INCLUDE 'GLK.h'
02555       INTEGER           kax,kay,nchx
02556       DOUBLE PRECISION  yl,yu,yy(*)
02557 *
02558       INTEGER           IX(m_MaxNb),IY(m_MaxNb)
02559       DOUBLE PRECISION  ai0,dx,aj0,ds,facy,aj,ai
02560       INTEGER           ipnt,i,iter,ipoin,irad1,irad2
02561       DOUBLE PRECISION  GLK_AproF
02562 *---------------------------------------------------
02563 
02564 * ...various types of line
02565 * ...distance between points is DS, radius of a point is IRAD
02566       IRAD2=6
02567       IRAD1=3
02568 * .............
02569       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02570      $  m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
02571       WRITE(m_ltx,'(A)') '% ========== plotting primitives =========='
02572       IF(m_tline .EQ. 1) THEN
02573 *   small filled circle
02574        DS = 10
02575        WRITE(m_ltx,'(8A,I3,A)')
02576      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02577      $   m_BS,'put(#1,#2){',m_BS,'circle*{',IRAD1,'}}}'
02578       ELSEIF(m_tline .EQ. 2) THEN
02579 *   small open circle
02580        DS = 10
02581        WRITE(m_ltx,'(8A,I3,A)')
02582      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02583      $   m_BS,'put(#1,#2){',m_BS,'circle{',IRAD1,'}}}'
02584       ELSEIF(m_tline .EQ. 3) THEN
02585 *   big filled circle
02586        DS = 20
02587        WRITE(m_ltx,'(8A,I3,A)')
02588      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02589      $   m_BS,'put(#1,#2){',m_BS,'circle*{',IRAD2,'}}}'
02590       ELSEIF(m_tline .EQ. 4) THEN
02591 *   big open circle
02592        DS = 20
02593        WRITE(m_ltx,'(8A,I3,A)')
02594      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02595      $   m_BS,'put(#1,#2){',m_BS,'circle{',IRAD2,'}}}'
02596 * Other symbols
02597       ELSEIF(m_tline .EQ. 5) THEN
02598        DS = 20
02599        WRITE(m_ltx,'(10A)')
02600      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02601      $   m_BS,'put(#1,#2){',m_BS,'makebox(0,0){$',m_BS,'diamond$}}}'
02602       ELSE
02603        DS = 20
02604        WRITE(m_ltx,'(10A)')
02605      $   m_BS,'newcommand{',m_BS,'R}[2]{',
02606      $   m_BS,'put(#1,#2){',m_BS,'makebox(0,0){$',m_BS,'star$}}}'
02607       ENDIF
02608       FACY = kay/(YU-YL)
02609 * plot first point
02610       AI  = 0.
02611       AJ  = (GLK_AproF( (AI/kax)*NCHX+0.5d0, NCHX, YY) -YL)*FACY
02612       IPNT =1
02613       IX(IPNT) = INT(AI)
02614       IY(IPNT) = INT(AJ)
02615       DX =  DS
02616       AI0 = AI
02617       AJ0 = AJ
02618 * plot next points
02619       DO 100 IPOIN=2,3000
02620 * iteration to get (approximately) equal distance among ploted points
02621       DO  50 ITER=1,3
02622       AI  = AI0+DX
02623       AJ  = (GLK_AproF( (AI/kax)*NCHX+0.5d0, NCHX, YY) -YL)*FACY
02624       DX  = DX *DS/SQRT(DX**2 + (AJ-AJ0)**2)
02625   50  CONTINUE
02626       IF(INT(AJ) .GE. 0.AND.INT(AJ) .LE. kay.AND.INT(AI) .LE. kax) THEN
02627          IPNT = IPNT+1
02628          IX(IPNT) = INT(AI)
02629          IY(IPNT) = INT(AJ)
02630       ENDIF
02631       AI0 = AI
02632       AJ0 = AJ
02633       IF(INT(AI) .GT. kax) GOTO 101
02634  100  CONTINUE
02635  101  CONTINUE
02636       WRITE(m_ltx,7000) (m_BS,IX(I),IY(I), I=1,IPNT)
02637 7000  FORMAT(4(A1,2HR{,I4,2H}{,I4,1H}:1X ))
02638       WRITE(m_ltx,'(2A)') m_BS,'end{picture}} % end of plotting line'
02639 * change line-style
02640       m_tline= m_tline+1
02641       IF(m_tline .GT. 2) m_tline=1
02642       END
02643 
02644       DOUBLE PRECISION  FUNCTION GLK_AproF(px,nch,yy)
02645 *     ************************************************
02646 * PX is a continuous extension of the m_index in array YY
02647       IMPLICIT NONE
02648       INTEGER           nch,ip
02649       DOUBLE PRECISION  px,yy(*),X,p
02650 *-----------------------------------------------------
02651       X=PX
02652       IF(X .LT. 0.0.OR.X .GT. FLOAT(NCH+1)) THEN
02653         GLK_AproF= -1E-20
02654         RETURN
02655       ENDIF
02656       IP=INT(X)
02657       IF(IP .LT. 2)     IP=2
02658       IF(IP .GT. NCH-2) IP=NCH-2
02659       P=X-IP
02660       GLK_AproF =
02661      $     -(1./6.)*P*(P-1)*(P-2)  *YY(IP-1)
02662      $     +(1./2.)*(P*P-1)*(P-2)  *YY(IP  )
02663      $     -(1./2.)*P*(P+1)*(P-2)  *YY(IP+1)
02664      $     +(1./6.)*P*(P*P-1)      *YY(IP+2)
02665       END
02666 
02667       SUBROUTINE GLK_PlTitle(title)
02668 *     *****************************
02669       IMPLICIT NONE
02670       INCLUDE 'GLK.h'
02671       SAVE
02672       CHARACTER*80 title
02673 *----------------------------------------
02674       m_KeyTit=1
02675       CALL GLK_Copch(title,m_titch(1))
02676       END
02677 
02678       SUBROUTINE GLK_PlCapt(lines)
02679 *     ****************************
02680 * This routine defines caption and should be called
02681 * before CALL GLK_Plot2, GLK_PlTable or bpltab2
02682 * The matrix CHARACTER*80 lines containes text of the caption ended
02683 * with the last line '% end-of-caption'
02684       IMPLICIT NONE
02685       CHARACTER*80 lines(*)
02686       INCLUDE 'GLK.h'
02687       SAVE
02688       INTEGER i
02689 *----------------------------------
02690       m_KeyTit=0
02691       DO i=1,m_titlen
02692          m_titch(i)=lines(i)
02693          m_KeyTit= m_KeyTit+1
02694          IF(lines(i) .EQ. '% end-of-caption' ) GOTO 100
02695       ENDDO
02696       CALL GLK_Retu1(' WARNING from GLK_PlCapt: to many lines =',m_titlen)
02697  100  CONTINUE
02698       END
02699 
02700       SUBROUTINE GLK_PlLabel(lines)
02701 *     *****************************
02702 * This should be envoked after 'CALL GLK_Plot2'
02703 * to add lines of TeX to a given plot
02704 *     ***********************************
02705       IMPLICIT NONE
02706       CHARACTER*80 lines(*)
02707       INCLUDE 'GLK.h'
02708       SAVE
02709       INTEGER i
02710 *----------------------------------
02711       m_KeyTit=0
02712       DO i=1,m_titlen
02713          m_titch(i)=lines(i)
02714          m_KeyTit= m_KeyTit+1
02715          IF(lines(i) .EQ. '% end-of-label' ) GOTO 100
02716       ENDDO
02717       CALL GLK_Retu1(' WARNING from GLK_PlLabel: to many lines =',m_titlen)
02718  100  CONTINUE
02719 *------------------------------!
02720 *   erase Ending               !
02721 *------------------------------!
02722       BACKSPACE(m_ltx)
02723       BACKSPACE(m_ltx)
02724 *
02725       DO i=1,m_KeyTit
02726         WRITE(m_ltx,'(A)')     m_titch(i)
02727       ENDDO
02728 *------------------------------!
02729 *   restore Ending             !
02730 *------------------------------!
02731       WRITE(m_ltx,'(2A)') m_BS,'end{picture} % close entire picture '
02732       IF(ABS(m_lint) .EQ. 2) THEN
02733          WRITE(m_ltx,'(A)') '%====== end of GLK_PlLabel =========='
02734       ELSE
02735          WRITE(m_ltx,'(2A)') m_BS,'end{figure}'
02736       ENDIF
02737       END
02738 
02739 
02740       SUBROUTINE GLK_Plot2(id,ch1,ch2,chmark,chxfmt,chyfmt)
02741 *     *****************************************************
02742 * The new, more user-friendly, version of older GLK_Plot
02743 * INPUT:
02744 *    ID          histogram identifier
02745 *    ch1 = ' '   normal new plot
02746 *        = 'S'   impose new plot on previous one
02747 *    ch2 = ' '   ploting line default, contour
02748 *        = '*'   error bars in midle of the bin
02749 *        = 'R'   error bars at Right edge of the bin
02750 *        = 'L'   error bars at Left  edge of the bin
02751 *        = 'C'   slanted continuous smooth line
02752 *    chmark =    TeX symbol for ploting points
02753 *    chxfmt =    format (string) for labeling x-axis
02754 *    chyfmt =    format (string) for labeling y-axis
02755 * Furthermore:
02756 * Captions are defined by means of
02757 *    CALL GLK_PlCapt(capt) before CALL GLK_Plot2
02758 *    where CHARACTER*80 capt(50) is content of
02759 *    caption, line by line, see also comments in GLK_PlCapt routine.
02760 * Additional text as a TeX source text can be appended by means of
02761 *    CALL GLK_PlLabel(lines) after CALL GLK_Plot2
02762 *    where CHARACTER*80 lines(50) is the TeX add-on.
02763 *    This is to be used to decorate plot with
02764 *    any kind marks, special labels and text on the plot.
02765 *
02766 *     ************************************
02767       IMPLICIT NONE
02768       INTEGER id
02769       CHARACTER ch1,ch2,chmark*(*)
02770       CHARACTER*8 chxfmt,chyfmt
02771       INCLUDE 'GLK.h'
02772       SAVE
02773       DOUBLE PRECISION  yy(m_MaxNb),yer(m_MaxNb)
02774       CHARACTER*80 title
02775 *---------------------------------------------------------------------
02776       LOGICAL GLK_Exist
02777       INTEGER kax,kay,incr,ker,nchx
02778       INTEGER iopsla,ioplog,ioperb,iopsc1,iopsc2,idum
02779       DOUBLE PRECISION   dxl,dxu,xu,xl,yu,yl
02780       CHARACTER chr
02781       DATA CHR /' '/
02782 * TeX Names of the error-bar command and of the point-mark command
02783       CHARACTER*1 chre, chrp1
02784       PARAMETER ( chre = 'E', chrp1= 'R' )
02785       CHARACTER*2 chrp
02786 * TeX Name of the point-mark command
02787       CHARACTER*1 chrx(12)
02788       DATA  chrx /'a','b','c','d','f','g','h','i','j','k','l','m'/
02789 *---------------------------------------------------------------------
02790 * RETURN if histo non-existing
02791       IF(.NOT.GLK_Exist(id)) GOTO 900
02792 * ...unpack histogram
02793       CALL GLK_UnPak(id,yy ,'    ',idum)
02794       CALL GLK_UnPak(id,yer,'ERRO',idum)
02795       CALL GLK_hinbo1(id,title,nchx,dxl,dxu)
02796 * Header
02797       kax=1200
02798       kay=1200
02799       IF(CH1 .EQ. 'S') THEN
02800 * Superimpose plot
02801         incr=incr+1
02802         BACKSPACE(m_ltx)
02803         BACKSPACE(m_ltx)
02804       ELSE
02805 * New frame only
02806         incr=1
02807         CHR=CH1
02808         CALL GLK_PlFrame(id,kax,kay,chxfmt,chyfmt)
02809 * The Y-range from first plot is preserved
02810         CALL GLK_Range1(id,yl,yu)
02811       ENDIF
02812 * The X-range as in histo
02813       xl = dxl
02814       xu = dxu
02815 *
02816       chrp= chrp1//chrx(incr)
02817       WRITE(m_ltx,'(A)')    '%=GLK_Plot2:  next plot (line) =========='
02818       WRITE(m_ltx,'(A,I10)')'%====HISTOGRAM ID=',ID
02819       WRITE(m_ltx,'(A,A70 )') '% ',TITLE
02820       CALL GLK_OptOut(id,ioplog,iopsla,ioperb,iopsc1,iopsc2)
02821       ker = ioperb-1
02822 * Default line type
02823       IF (iopsla .EQ. 2) THEN
02824          CHR='C'
02825       ELSE
02826          CHR=' '
02827       ENDIF
02828 * User defined line-type
02829       IF (CH2 .EQ. 'B')   CHR=' '
02830 *...marks in the midle of the bin
02831       IF (CH2 .EQ. '*')   CHR='*'
02832 *...marks on the right edge of the bin
02833       IF (CH2 .EQ. 'R')   CHR='R'
02834 *...marks on the left edge of the bin
02835       IF (CH2 .EQ. 'L')   CHR='L'
02836       IF (CH2 .EQ. 'C')   CHR='C'
02837 *...various types of lines
02838       IF     (CHR .EQ. ' ') THEN
02839 *...contour line used for histogram
02840           CALL GLK_PlKont(kax,kay,nchx,yl,yu,yy,ker,yer)
02841       ELSE IF(CHR .EQ. '*' .OR. CHR .EQ. 'R'.OR. CHR .EQ. 'L') THEN
02842 *...marks on the right/left/midle of the bin
02843          CALL GLK_PlMark(kax,kay,nchx,yl,yu,yy,ker,yer,chmark,chr,chrp,chre)
02844       ELSE IF(CHR .EQ. 'C') THEN
02845 *...slanted (dotted) line in plotting non-MC functions
02846           CALL GLK_PlCirc(kax,kay,nchx,yl,yu,yy)
02847       ENDIF
02848 *------------------------------!
02849 *        ENDing                !
02850 *------------------------------!
02851       WRITE(m_ltx,'(2A)') m_BS,'end{picture} % close entire picture '
02852       IF(ABS(m_lint) .EQ. 2) THEN
02853          WRITE(m_ltx,'(A)') '%== GLK_Plot2:  end of plot  =========='
02854       ELSE
02855          WRITE(m_ltx,'(2A)') m_BS,'end{figure}'
02856       ENDIF
02857       RETURN
02858   900 CALL GLK_Stop1('+++GLK_Plot2: Nonexistig histo, skipped, id= ',ID)
02859       END
02860 
02861       SUBROUTINE GLK_PlFrame(id,kax,kay,chxfmt,chyfmt)
02862 *     ************************************************
02863       IMPLICIT NONE
02864       INTEGER id,kax,kay
02865       CHARACTER chxfmt*(*),chyfmt*(*)
02866       INCLUDE 'GLK.h'
02867       SAVE
02868 *---------------------------------------------------
02869       CHARACTER*80 title
02870       DOUBLE PRECISION    dxl,dxu,xl,xu,yl,yu
02871       INTEGER  icont,i,nchx
02872       DATA icont/0/
02873 *---------------------------------------------------
02874       icont=icont+1
02875       CALL GLK_hinbo1(id,title,nchx,dxl,dxu)
02876       xl = dxl
02877       xu = dxu
02878       CALL GLK_Range1(id,yl,yu)
02879 *
02880       IF(icont .GT. 1) WRITE(m_ltx,'(2A)') m_BS,'newpage'
02881 *------------------------------!
02882 *           Header
02883 *------------------------------!
02884       WRITE(m_ltx,'(A)') ' '
02885       WRITE(m_ltx,'(A)') ' '
02886       WRITE(m_ltx,'(A)')
02887      $'%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%'
02888       WRITE(m_ltx,'(A)')
02889      $'%%%%%%%%%%%%%%%%%%%%%%GLK_PlFrame%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%'
02890       IF(ABS(m_lint) .EQ. 2) THEN
02891          WRITE(m_ltx,'(2A)') m_BS,'noindent'
02892       ELSE
02893          WRITE(m_ltx,'(2A)') m_BS,'begin{figure}[!ht]'
02894          WRITE(m_ltx,'(2A)') m_BS,'centering'
02895          WRITE(m_ltx,'(2A)') m_BS,'htmlimage{scale=1.4}'
02896       ENDIF
02897 *------------------------------!
02898 * General Caption
02899 *------------------------------!
02900       IF(ABS(m_lint) .NE. 2) THEN
02901          WRITE(m_ltx,'(6A)')
02902      $        m_BS,'caption{',m_BS,'footnotesize',m_BS,'sf'
02903          DO i=1,m_KeyTit
02904             WRITE(m_ltx,'(A)')     m_titch(i)
02905          ENDDO
02906          WRITE(m_ltx,'(A)') '}'
02907       ENDIF
02908 *------------------------------!
02909 * Frames and labels
02910 *------------------------------!
02911       WRITE(m_ltx,'(A)') '% =========== big frame, title etc. ======='
02912       WRITE(m_ltx,'(4A)') m_BS,'setlength{',m_BS,'unitlength}{0.1mm}'
02913       WRITE(m_ltx,'(2A)') m_BS,'begin{picture}(1600,1500)'
02914       IF( m_lint .LT. 0) THEN
02915 * Big frame usefull for debuging
02916          WRITE(m_ltx,'(4A)')
02917      $        m_BS,'put(0,0){',m_BS,'framebox(1600,1500){ }}'
02918       ENDIF
02919       WRITE(m_ltx,'(A)') '% =========== small frame, labeled axis ==='
02920       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02921      $    m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
02922       WRITE(m_ltx,'(4A,I4,A,I4,A)')
02923      $    m_BS,'put(0,0){',m_BS,'framebox( ',kax,',',kay,'){ }}'
02924       WRITE(m_ltx,'(A)') '% =========== x and y axis ================'
02925       CALL GLK_AxisX(kax,xl,xu,chxfmt)
02926       CALL GLK_AxisY(kay,yl,yu,chyfmt)
02927       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}'
02928      $                ,'% end of plotting labeled axis'
02929       END
02930 
02931       SUBROUTINE GLK_AxisX(kay,yl,yu,chxfmt)
02932 *     ***************************************
02933 * plotting x-axis with long and short tips
02934       IMPLICIT NONE
02935       INTEGER  kay
02936       DOUBLE PRECISION    yl,yu
02937       CHARACTER chxfmt*16
02938       INCLUDE 'GLK.h'
02939       SAVE
02940 *-------------------------------------------------------
02941       CHARACTER*64 fmt1,fmt2
02942       PARAMETER (fmt1 = '(2A,F8.2,A,F8.2,A,I4,3A)')
02943       PARAMETER (fmt2 = '(2A,F8.2,A,I4,A,F8.2,A,I4,3A)')
02944       DOUBLE PRECISION   dy,ddy,ddyl,yy0l,ddys,yy0s,p0s,pds,scmx,p0l,pdl
02945       INTEGER ly,jy,nlt,nts,lex,k,n
02946       DOUBLE PRECISION  tipsy(20)
02947 *-------------------------------------------------------
02948       dy= ABS(yu-yl)
02949       ly = NINT( LOG10(dy) -0.4999999d0 )
02950       jy = NINT(dy/10d0**ly)
02951       ddyl = dy*10d0**(-ly)
02952       IF( jy .EQ. 1)                 ddyl = 10d0**ly*0.25d0
02953       IF( jy .GE. 2 .AND. jy .LE. 3) ddyl = 10d0**ly*0.5d0
02954       IF( jy .GE. 4 .AND. jy .LE. 6) ddyl = 10d0**ly*1.0d0
02955       IF( jy .GE. 7)                 ddyl = 10d0**ly*2.0d0
02956       WRITE(m_ltx,'(A)') '% -------GLK_AxisX---- '
02957       WRITE(m_ltx,'(A,I4)') '%  JY= ',JY
02958 *-------
02959       nlt = INT(dy/ddyl)
02960       nlt = MAX0(MIN0(nlt,20),1)+1
02961       yy0l = NINT(yl/ddyl+0.5d0)*ddyl
02962       ddys = ddyl/10d0
02963       yy0s = NINT(yl/ddys+0.4999999d0)*ddys
02964       p0l = kay*(yy0l-yl)/(yu-yl)
02965       pdl = kay*ddyl/(yu-yl)
02966       p0s = kay*(yy0s-yl)/(yu-yl)
02967       pds = kay*ddys/(yu-yl)
02968       nlt = INT(ABS(yu-yy0l)/ddyl+0.0000001d0)+1
02969       nts = INT(abs(yu-yy0s)/ddys+0.0000001d0)+1
02970       DO n=1,nlt
02971          tipsy(n) =yy0l+ ddyl*(n-1)
02972       ENDDO
02973       WRITE(m_ltx,fmt1)
02974      $ m_BS,'multiput('  ,P0L,  ',0)('  ,PDL,  ',0){'  ,NLT,  '}{',
02975      $ m_BS,'line(0,1){25}}',
02976      $ m_BS,'multiput('  ,P0S,  ',0)('  ,PDS,  ',0){'  ,NTS,  '}{',
02977      $ m_BS,'line(0,1){10}}'
02978       WRITE(m_ltx,fmt2)
02979      $ m_BS,'multiput('  ,P0L,  ','  ,kay,  ')('  ,PDL,  ',0){'  ,NLT,
02980      $ '}{'  ,m_BS,  'line(0,-1){25}}',
02981      $ m_BS,'multiput('  ,P0S,  ','  ,kay,  ')('  ,PDS,  ',0){'  ,NTS,
02982      $ '}{'  ,m_BS,  'line(0,-1){10}}'
02983 * ...labeling of axis
02984       scmx = DMAX1(DABS(yl),DABS(YU))
02985       lex  = NINT( LOG10(scmx) -0.50001)
02986       DO n=1,nlt
02987          k = nint(kay*(tipsy(n)-yl)/(yu-yl))
02988          IF(lex .LE. 3 .AND. lex .GE. -3) THEN
02989 * ...without exponent
02990            WRITE(m_ltx,'(2A,I4,5A,'//chxfmt//',A)')
02991      $     m_BS,'put(',K,',-25){',m_BS,'makebox(0,0)[t]{',
02992      $           m_BS,'Large $ ', TIPSY(N), ' $}}'
02993          ELSE
02994 * ...with exponent
02995            WRITE(m_ltx,'(2A,I4,5A,'//chxfmt//',2A,I4,A)')
02996      $     m_BS,'put('  ,K,  ',-25){',m_BS,'makebox(0,0)[t]{',
02997      $     m_BS,'Large $ ',
02998      $     TIPSY(N)/(10d0**LEX),m_BS,'cdot 10^{',LEX,'} $}}'
02999          ENDIF
03000       ENDDO
03001       END
03002 
03003       SUBROUTINE GLK_AxisY(kay,yl,yu,chyfmt)
03004 *     ***************************************
03005 * plotting y-axis with long and short tips
03006       IMPLICIT NONE
03007       INTEGER  kay
03008       DOUBLE PRECISION    yl,yu
03009       CHARACTER chyfmt*16
03010       INCLUDE 'GLK.h'
03011       SAVE
03012       DOUBLE PRECISION  tipsy(20)
03013 *------------------------------------------------------------------
03014       CHARACTER*64 fmt1,fmt2
03015       PARAMETER (fmt1 = '(2A,F8.2,A,F8.2,A,I4,3A)')
03016       PARAMETER (fmt2 = '(2A,I4,A,F8.2,A,F8.2,A,I4,3A)')
03017       INTEGER ly,jy,nlt,nts,lex,n,k
03018       DOUBLE PRECISION   ddyl,dy,yy0l,p0l,pdl,pds,scmx,z0l,p0s,yy0s,ddys
03019 *------------------------------------------------------------------
03020       dy= ABS(yu-yl)
03021       ly = NINT( log10(dy) -0.49999999d0 )
03022       jy = NINT(dy/10d0**ly)
03023       ddyl = dy*10d0**(-ly)
03024       IF( jy .EQ. 1)                 ddyl = 10d0**ly*0.25d0
03025       IF( jy .GE. 2 .AND. jy .LE. 3) ddyl = 10d0**ly*0.5d0
03026       IF( jy .GE. 4 .AND. jy .LE. 6) ddyl = 10d0**ly*1.0d0
03027       IF( jy .GE. 7)                 ddyl = 10d0**ly*2.0d0
03028       WRITE(m_ltx,'(A)') '% --------GLK_SAxisY------- '
03029       WRITE(m_ltx,'(A,I4)') '%  JY= ',JY
03030 *-------
03031       nlt = INT(dy/ddyl)
03032       nlt = MAX0(MIN0(nlt,20),1)+1
03033       yy0l = NINT(yl/ddyl+0.4999999d0)*ddyl
03034       ddys = ddyl/10d0
03035       yy0s = nint(yl/ddys+0.5d0)*ddys
03036       p0l = kay*(yy0l-yl)/(yu-yl)
03037       pdl = kay*ddyl/(yu-yl)
03038       p0s = kay*(yy0s-yl)/(yu-yl)
03039       pds = kay*ddys/(yu-yl)
03040       nlt= INT(ABS(yu-yy0l)/ddyl+0.0000001d0) +1
03041       nts= INT(ABS(yu-yy0s)/ddys+0.0000001d0) +1
03042       DO N=1,NLT
03043          tipsy(n) =yy0l+ ddyl*(n-1)
03044       ENDDO
03045 * plotting tics on vertical axis
03046       WRITE(m_ltx,fmt1)
03047      $ m_BS,'multiput(0,'  ,P0L,  ')(0,'  ,PDL  ,'){'  ,NLT,  '}{', m_BS,'line(1,0){25}}',
03048      $ m_BS,'multiput(0,'  ,P0S,  ')(0,'  ,PDS,  '){'  ,NTS,  '}{', m_BS,'line(1,0){10}}'
03049       WRITE(m_ltx,fmt2)
03050      $ m_BS,'multiput('  ,kay,  ','  ,P0L,  ')(0,'  ,PDL,  '){'  ,NLT,
03051      $ '}{',m_BS,'line(-1,0){25}}',
03052      $ m_BS,'multiput('  ,kay,  ','  ,P0S,  ')(0,'  ,PDS,  '){'  ,NTS,
03053      $ '}{',m_BS,'line(-1,0){10}}'
03054 * ...Zero line if necessary
03055       Z0L = kay*(-YL)/(YU-YL)
03056       IF( (Z0L .GT. 0D0) .AND. (Z0L .LT. FLOAT(kay)) )
03057      $ WRITE(m_ltx,'(2A,F8.2,3A,I4,A)') m_BS,'put(0,'  ,Z0L,  '){',m_BS,'line(1,0){'  ,kay,  '}}'
03058 * ...labeling of axis
03059       SCMX = DMAX1(DABS(YL),DABS(YU))
03060       LEX  = NINT( LOG10(SCMX) -0.50001d0)
03061       DO n=1,nlt
03062          k = nint(kay*(tipsy(n)-yl)/(yu-yl))
03063          IF(lex .LE. 3 .AND. lex .GE. -3) THEN
03064 * ...without exponent
03065             WRITE(m_ltx,'(2A,I4,5A,'//chyfmt//',A)')
03066      $           m_BS,'put(-25,'  ,K,  '){',m_BS,'makebox(0,0)[r]{',
03067      $           m_BS,'Large $ '  ,TIPSY(N),  ' $}}'
03068          ELSE
03069 * ...with exponent
03070             WRITE(m_ltx,'(2A,I4,5A,'//chyfmt//',2A,I4,A)')
03071      $           m_BS,'put(-25,'  ,K,  '){',m_BS,'makebox(0,0)[r]{',
03072      $           m_BS,'Large $ ',
03073      $           TIPSY(N)/(10d0**LEX),  m_BS,'cdot 10^{'  ,LEX,  '} $}}'
03074       ENDIF
03075       ENDDO
03076       END
03077 
03078       SUBROUTINE GLK_PlKont(kax,kay,nchx,yl,yu,yy,ker,yer)
03079 */////////////////////////////////////////////////////////////////////////////////////
03080 *//                                                                                 //
03081 *//             Plotting contour line for histogram (formely PlHis)                 //
03082 *//                                                                                 //
03083 */////////////////////////////////////////////////////////////////////////////////////
03084       IMPLICIT NONE
03085       INTEGER kax,kay,nchx,ker
03086       DOUBLE PRECISION   yl, yu, yy(*),yer(*),z0l
03087       INCLUDE 'GLK.h'
03088       SAVE
03089 *---------------------------------------------------
03090       CHARACTER*80 fmt1
03091       INTEGER  ix0,iy0,ib,ix1,iy1,ie,ierr,ix2,idy,idx
03092       DOUBLE PRECISION    yib
03093 *---------------------------------------------------
03094       WRITE(m_ltx,'(4A,I4,A,I4,A)') m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
03095       WRITE(m_ltx,'(A)') '% ========== plotting primitives =========='
03096 * Color string, optionaly
03097       IF(m_KeyCol .EQ. 1) THEN
03098          WRITE(m_ltx,'(A)') m_Color
03099          m_KeyCol = 0
03100       ENDIF
03101 *...short macros for vertical/horizontal straight lines
03102       WRITE(m_ltx,'(8A)')
03103      $     m_BS,'newcommand{',m_BS,'x}[3]{',m_BS,'put(#1,#2){', m_BS,'line(1,0){#3}}}'
03104       WRITE(m_ltx,'(8A)')
03105      $     m_BS,'newcommand{',m_BS,'y}[3]{',m_BS,'put(#1,#2){', m_BS,'line(0,1){#3}}}'
03106       WRITE(m_ltx,'(8A)')
03107      $     m_BS,'newcommand{',m_BS,'z}[3]{',m_BS,'put(#1,#2){', m_BS,'line(0,-1){#3}}}'
03108 *   error bars
03109       WRITE(m_ltx,'(8A)')
03110      $     m_BS,'newcommand{',m_BS,'e}[3]{', m_BS,'put(#1,#2){',m_BS,'line(0,1){#3}}}'
03111 * Starting point for the line
03112       ix0=0
03113       iy0=0
03114 * Start at Zero line if possible
03115       z0l = kay*(-yl)/(yu-yl)
03116       IF( (z0l .GT. 0d0) .AND. (z0l .LT. FLOAT(kay)) )  iy0=z0l
03117       DO ib=1,nchx
03118          yib = yy(ib)
03119          ix1 = NINT(kax*(ib-0.00001d0)/nchx) ! new x
03120          iy1 = NINT(kay*(yib-yl)/(yu-yl))    ! new y
03121          iy1 = MIN(MAX(iy1,-1),kay+1)        ! cosmetics
03122          idx = ix1-ix0                       ! delta x
03123          idy = iy1-iy0                       ! delta y
03124          fmt1 = '(2(2a,i4,a,i4,a,i4,a))'
03125          IF(iy1 .GE. 0 .AND. iy1 .LE. kay) THEN
03126             IF( idy .GE. 0) THEN             ! up
03127                WRITE(m_ltx,fmt1) m_BS,'y{',ix0,'}{',iy0,'}{',idy,'}',
03128      $                           m_BS,'x{',ix0,'}{',iy1,'}{',idx,'}'
03129             ELSE                             ! down
03130                WRITE(m_ltx,fmt1) m_BS,'z{',IX0,'}{',IY0,'}{',-idy,'}',
03131      $                           m_BS,'x{',IX0,'}{',IY1,'}{',idx,'}'
03132             ENDIF
03133          ENDIF
03134          ix0=ix1
03135          iy0=iy1
03136          IF(ker .EQ. 1) THEN
03137             ix2  = NINT(kax*(ib-0.5000d0)/nchx)
03138             ierr = NINT(kay*((yy(ib)-yer(ib))-yl)/(yu-yl))  ! bottom of error bar
03139             ie = NINT(kay*yer(ib)/(yu-yl))                  ! total length of error bar
03140 *        Cosmetics
03141             IF(ierr .LT. 0) THEN
03142                ie= ie+ierr
03143                ierr = 0
03144             ENDIF
03145             IF( (ierr+2*ie) .GT. kay) THEN
03146                ie= IABS(kay-ierr)/2
03147             ENDIF
03148             IF( (iy1.GE.0).AND.(iy1.LE. kay).AND.(ABS(1d0*ierr).LE.9999d0).AND.(2d0*ie.LE.9999d0) ) 
03149      $           WRITE(m_ltx,8000) m_BS,ix2,ierr,2*ie
03150          ENDIF
03151       ENDDO
03152 8000  FORMAT(4(A1,2He{,I4,2H}{,I5,2H}{,I4,1H}:1X ))
03153       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}', ' % end of plotting histogram'
03154 * change line-style
03155       m_tline= m_tline+1
03156       IF(m_tline .GT. 2) m_tline=1
03157       END
03158 
03159       SUBROUTINE GLK_PlMark(kax,kay,nchx,yl,yu,yy,ker,yer,chmark,chr,chr2,chr3)
03160 */////////////////////////////////////////////////////////////////////////////////////
03161 *//                                                                                 //
03162 *//                       marks in the midle of the bin                             //
03163 *//                                                                                 //
03164 */////////////////////////////////////////////////////////////////////////////////////
03165       IMPLICIT NONE
03166       INTEGER     kax,kay,nchx,ker
03167       DOUBLE PRECISION       yl,yu, yy(*),yer(*)
03168       CHARACTER*1 chr
03169       CHARACTER   chmark*(*),chr2*(*),chr3*(*)
03170 *---------------------------------------------------
03171       INCLUDE 'GLK.h'
03172       SAVE
03173       INTEGER    ib,ix1,iy1,ierr,ie
03174 *---------------------------------------------------
03175       WRITE(m_ltx,'(4A,I4,A,I4,A)') m_BS,'put(300,250){',m_BS,'begin{picture}( ',kax,',',kay,')'
03176       WRITE(m_ltx,'(A)') '% ===GLK_PlMark: plotting primitives ======'
03177 * Color string, optionaly
03178       IF(m_KeyCol .EQ. 1) THEN
03179          WRITE(m_ltx,'(A)') m_Color
03180          m_KeyCol = 0
03181       ENDIF
03182 * Plotting symbol
03183       WRITE(m_ltx,'(10A)') m_BS,'newcommand{',m_BS,chr2  , '}[2]{', m_BS,'put(#1,#2){',chmark,'}}'
03184 * Error bar symbol
03185       WRITE(m_ltx,'(10A)')
03186      $   m_BS,'newcommand{',m_BS,chr3  , '}[3]{', m_BS,'put(#1,#2){',m_BS,'line(0,1){#3}}}'
03187 
03188       DO ib=1,nchx
03189          IF(chr .EQ. '*') THEN
03190             ix1 = NINT(kax*(ib-0.5000d0)/nchx) ! Midle of bin
03191          ELSEIF(chr .EQ. 'R') THEN
03192             ix1 = NINT(kax*(ib*1d0)/nchx)      ! Right edge of bin
03193          ELSEIF(chr .EQ. 'L') THEN
03194             ix1 = NINT(kax*(ib-1d0)/nchx)      ! Left edge of bin
03195          ELSE
03196             WRITE(6,*) '+++++ plamark: wrong line type:',chr
03197             RETURN
03198          ENDIF
03199          iy1 = NINT(kay*(yy(ib)-yl)/(yu-yl))
03200          IF(iy1 .GE. 0 .AND. iy1 .LE. kay)
03201      $   WRITE(m_ltx,'(A,A,A,I4,A,I4,A)')
03202      $               m_BS,chr2, '{' ,IX1, '}{' ,IY1, '}'
03203          IF(ker .EQ. 1) THEN
03204             ierr = NINT(kay*((yy(ib)-yer(ib))-yl)/(yu-yl)) ! bottom of error bar
03205             ie   = NINT(kay*yer(ib)/(yu-yl))               ! total length of error bar
03206 *        Cosmetics
03207             IF(ierr .LT. 0) THEN
03208                ie= ie+ierr
03209                ierr = 0
03210             ENDIF
03211             IF( (ierr+2*ie) .GT. kay) THEN
03212                ie= IABS(kay-ierr)/2
03213             ENDIF
03214             IF((iy1.GE.0) .AND.(iy1.LE.kay) .AND.(ABS(1d0*ierr).LE.9999d0) .AND.(2d0*ie.LE.9999d0))
03215      $      WRITE(m_ltx,'(A,A,A,I4,A,I5,A,I4,A)')
03216      $          m_BS, chr3,  '{'  ,IX1, '}{'  ,ierr, '}{'  ,2*ie,   '}'
03217          ENDIF
03218       ENDDO
03219       WRITE(m_ltx,'(3A)') m_BS,'end{picture}}',
03220      $    ' % end of plotting histogram'
03221       END
03222 
03223 
03224       SUBROUTINE GLK_PlTable(Npl,idl,capt,fmt,nch1,incr,npag)
03225 *     ******************************************************
03226 * Tables in TeX, up to 9 columns
03227 * Npl           = numbers of columns/histograms
03228 * idl(1:Npl)    = list of histo id's
03229 * capt(1:Npl+1) = list of captions above each column
03230 * fmt(1:1)      = format to print x(i) in first columb,
03231 *                 h(i) and error he(i) in further columns
03232 * nch1,incr     = raws are printet in the sequence
03233 *                 (h(i),he(i),i=nch1,nbin,incr), nbin is no. of bins.
03234 * npag          = 0 no page eject, =1 with page eject
03235 *     ******************************************************
03236       IMPLICIT NONE
03237 *--------------- parameters ------------
03238       INTEGER        Npl,idl(*),nch1,incr,npag
03239       CHARACTER*(*)  capt(*)
03240       CHARACTER*(*)  fmt(3)
03241 *-------------------------------------------
03242       INCLUDE 'GLK.h'
03243       SAVE
03244 *---------------------------------------------------
03245       CHARACTER*16 fmt1,fmt2,fmt3
03246       LOGICAL GLK_Exist
03247       INTEGER   i,j,k,n,nchx,nplt,idum,id1,id
03248       INTEGER   iopsc1,ioperb,iopsla,iopsc2,ioplog
03249       DOUBLE PRECISION     xl,xu,dxl,dxu,xi
03250       DOUBLE PRECISION     yyy(m_MaxNb),yer(m_MaxNb),bi(m_MaxNb,9),er(m_MaxNb,9)
03251       CHARACTER*80 title
03252       CHARACTER*1 Cn(9)
03253       DATA Cn /'1','2','3','4','5','6','7','8','9'/
03254 *-----------------------------------------------------------------------------
03255 * Return if histo non-existing or to many columns
03256       IF(.NOT.GLK_EXIST(ID)) GOTO 900
03257       IF(Npl .GT. 9 )     GOTO 901
03258       fmt1 = fmt(1)
03259       fmt2 = fmt(2)
03260       fmt3 = fmt(3)
03261 *
03262 * npack histograms
03263       id1=idl(1)
03264       CALL GLK_hinbo1( id1,title,nchx,dxl,dxu)
03265       xl = dxl
03266       xu = dxu
03267       DO n=1,Npl
03268         CALL GLK_UnPak( idl(n),yyy ,'    ',idum)
03269         CALL GLK_UnPak( idl(n),yer ,'ERRO',idum)
03270         DO k=1,nchx
03271            bi(k,n)=yyy(k)
03272            er(k,n)=yer(k)
03273         ENDDO
03274       ENDDO
03275 *------------------------------!
03276 *           Header
03277 *------------------------------!
03278       WRITE(m_ltx,'(A)') ' '
03279       WRITE(m_ltx,'(A)') ' '
03280       WRITE(m_ltx,'(A)') '% ========================================='
03281       WRITE(m_ltx,'(A)') '% ============= begin table ==============='
03282       WRITE(m_ltx,'(2A)') m_BS,'begin{table}[!ht]'
03283       WRITE(m_ltx,'(2A)') m_BS,'centering'
03284 *------------------------------!
03285 * Central Caption
03286 *------------------------------!
03287       WRITE(m_ltx,'(4A)') m_BS,'caption{',m_BS,'small'
03288       DO i=1,m_KeyTit
03289         WRITE(m_ltx,'(A)')     m_titch(i)
03290       ENDDO
03291       WRITE(m_ltx,'(A)') '}'
03292 *------------------------------!
03293 * Tabular header
03294 *------------------------------!
03295       WRITE(m_ltx,'(20A)') m_BS,
03296 'begin{tabular}     $ {|',  ('|c',j=1,Npl+1),  '||}'
03297 *
03298       WRITE(m_ltx,'(4A)') m_BS,'hline',m_BS,'hline'
03299 *------------------------------!
03300 * Captions in columns
03301 *------------------------------!
03302       WRITE(m_ltx,'(2A)') capt(1),('&',capt(j+1),j=1,Npl)
03303 *
03304       WRITE(m_ltx,'(2A)') m_BS,m_BS
03305       WRITE(m_ltx,'(2A)') m_BS,'hline'
03306 *----------------------------------------!
03307 * Table content
03308 * Note that by default RIGHT EDGE of bin is printed, as necessary for
03309 * cumulative distributions, this can be changed with SLAN option
03310 *----------------------------------------!
03311       CALL GLK_OptOut(idl(1),ioplog,iopsla,ioperb,iopsc1,iopsc2)
03312       DO k=nch1,nchx,incr
03313         xi= dxl + (dxu-dxl)*k/(1d0*nchx)
03314         IF(iopsla.eq.2) xi= dxl + (dxu-dxl)*(k-0.5d0)/(1d0*nchx)
03315         IF(ioperb.eq.2) THEN
03316         WRITE(m_ltx,'(A,'//fmt1//','//Cn(Npl)//'(A,'//fmt2//',A,A,'//fmt3//'),  A)')
03317      $               '$', xi, ('$ & $', bi(k,j), m_BS, 'pm', er(k,j), j=1,Npl), '$'
03318         WRITE(m_ltx,'(2A)') m_BS,m_BS
03319         ELSE
03320         WRITE(m_ltx,'(A,'//fmt1//','//Cn(Npl)//'(A,'//fmt2//'),  A)')
03321      $               '$', xi, ('$ & $', bi(k,j), j=1,Npl), '$'
03322         WRITE(m_ltx,'(2A)') m_BS,m_BS
03323         ENDIF
03324       ENDDO
03325 *------------------------------!
03326 * Ending
03327 *------------------------------!
03328       WRITE(m_ltx,'(4A)') m_BS,'hline',m_BS,'hline'
03329       WRITE(m_ltx,'(2A)') m_BS,'end{tabular}'
03330       WRITE(m_ltx,'(2A)') m_BS,'end{table}'
03331       WRITE(m_ltx,'(A)') '% ============= end   table ==============='
03332       WRITE(m_ltx,'(A)') '% ========================================='
03333       IF(npag .NE. 0) WRITE(m_ltx,'(2A)') m_BS,'newpage'
03334 
03335       RETURN
03336  900  CALL GLK_Retu1('++++ GLK_PlTable: Nonexistig histo id=',ID)
03337       RETURN
03338  901  CALL GLK_Retu1('++++ GLK_PlTable: To many columns Nplt=',Nplt)
03339       END
03340 
03341       SUBROUTINE GLK_PlTable2(Npl,idl,ccapt,mcapt,fmt,chr1,chr2,chr3)
03342 *     ***************************************************************
03343 * Tables in TeX, up to 9 columns
03344 * Npl           = numbers of columns/histograms
03345 * idl(1:Npl)    = list of histo id's
03346 * ccapt(1:Npl+1)= list of column-captions above each column
03347 * mcapt         = multicolumn header, none if mcapt=' ',
03348 * fmt(1:1)      = format to print x(i) in first columb,
03349 *                 h(i) and error he(i) in further columns
03350 * chr1          = ' ' normal default, ='S' the Same table continued
03351 * chr2          = ' ' midle of the bin for x(i) in the first column
03352 *               = 'R' right edge,     ='L' left edge of the bin
03353 * chr3          = ' ' no page eject,  ='E' with page eject at the end.
03354 * Furthermore:
03355 * Captions are defined by means of
03356 *    CALL GLK_PlCapt(capt) before CALL GLK_PlTable2
03357 *    where CHARACTER*80 capt(50) is content of
03358 *    caption, line by line, see also comments in GLK_PlCapt routine.
03359 *
03360 *     ******************************************************
03361       IMPLICIT NONE
03362 *-------------- parameters--------------
03363       INTEGER       Npl,idl(*)
03364       CHARACTER*(*) ccapt(*)
03365       CHARACTER*(*) fmt(3)
03366       CHARACTER*1   chr1,chr2,chr3
03367       CHARACTER*(*) mcapt
03368 *----------------------------------------------------------------------
03369       INCLUDE 'GLK.h'
03370       SAVE
03371 *----------------------------------------------------------------------
03372       CHARACTER*16 fmt1,fmt2,fmt3
03373       LOGICAL GLK_Exist
03374       INTEGER   iopsc1,ioperb,iopsla,iopsc2,ioplog
03375       INTEGER   i,j,k,n,idum,id1,id,nchx,Nplt
03376       DOUBLE PRECISION     xl,xu,xi,dxu,dxl
03377       DOUBLE PRECISION     yyy(m_MaxNb),yer(m_MaxNb),bi(m_MaxNb,9),er(m_MaxNb,9)
03378       CHARACTER*80 title
03379       CHARACTER*1 Cn(9)
03380       INTEGER   k1,k2,k3
03381       DATA Cn /'1','2','3','4','5','6','7','8','9'/
03382 *----------------------------------------------------------------------
03383 * RETURN if histo non-existing or to many columns
03384       IF(.NOT.GLK_EXIST(ID)) GOTO 900
03385       IF(Npl .GT. 9 )     GOTO 901
03386       fmt1 = fmt(1)
03387       fmt2 = fmt(2)
03388       fmt3 = fmt(3)
03389 *
03390 * unpack histograms
03391       id1 = idl(1)
03392       CALL GLK_hinbo1( id1,title,nchx,dxl,dxu)
03393       xl = dxl
03394       xu = dxu
03395       DO n=1,Npl
03396          CALL GLK_UnPak( idl(n),yyy ,'    ',idum)
03397          CALL GLK_UnPak( idl(n),yer ,'ERRO',idum)
03398          DO k=1,nchx
03399             bi(k,n)=yyy(k)
03400             er(k,n)=yer(k)
03401          ENDDO
03402       ENDDO
03403 
03404       IF(chr1 .EQ. ' ' ) THEN
03405 *------------------------------!
03406 *           Header
03407 *------------------------------!
03408          WRITE(m_ltx,'(A)') ' '
03409          WRITE(m_ltx,'(A)') ' '
03410          WRITE(m_ltx,'(A)') '% ========================================'
03411          WRITE(m_ltx,'(A)') '% ============ begin table ==============='
03412 *
03413          IF(ABS(m_lint) .EQ. 2 ) THEN
03414             WRITE(m_ltx,'(2A)') m_BS,'noindent'
03415          ELSE
03416             WRITE(m_ltx,'(2A)') m_BS,'begin{table}[!ht]'
03417             WRITE(m_ltx,'(2A)') m_BS,'centering'
03418          ENDIF
03419 *------------------------------!
03420 * Central Caption
03421 *------------------------------!
03422          IF(ABS(m_lint) .NE. 2 ) THEN
03423             WRITE(m_ltx,'(6A)')
03424      $           m_BS,'caption{',m_BS,'footnotesize',m_BS,'sf'
03425             DO i=1,m_KeyTit
03426                WRITE(m_ltx,'(A)')     m_titch(i)
03427             ENDDO
03428             WRITE(m_ltx,'(A)') '}'
03429          ENDIF
03430 *------------------------------!
03431 * Tabular header
03432 *------------------------------!
03433          WRITE(m_ltx,'(20A)') m_BS,
03434 'begin{tabular}     $        {|',  ('|c',j=1,Npl+1),  '||}'
03435          WRITE(m_ltx,'(4A)') m_BS,'hline',m_BS,'hline'
03436 *------------------------------!
03437 * Captions in columns
03438 *------------------------------!
03439          WRITE(m_ltx,'(2A)') ccapt(1),('&',ccapt(j+1),j=1,Npl)
03440 *------------------------------!
03441 * Append previous table
03442 *------------------------------!
03443       ELSEIF(chr1 .EQ. 'S' ) THEN
03444          DO i=1,7
03445             BACKSPACE(m_ltx)
03446          ENDDO
03447       ELSE
03448          WRITE(*,*) ' ++++ GLK_PlTable2: WRONG chr1 ' ,chr1
03449       ENDIF
03450 
03451       WRITE(m_ltx,'(2A)') m_BS,m_BS
03452       WRITE(m_ltx,'(2A)') m_BS,'hline'
03453 
03454 *------------------------------!
03455 * Optional multicolumn caption
03456 *------------------------------!
03457       IF(mcapt .NE. ' ') THEN
03458          WRITE(m_ltx,'(3A,I2,A)') '& ',m_BS,'multicolumn{',Npl,'}{c||}{'
03459          WRITE(m_ltx,'(3A)') '     ',mcapt, ' }'
03460          WRITE(m_ltx,'(2A)') m_BS,m_BS
03461          WRITE(m_ltx,'(2A)') m_BS,'hline'
03462       ENDIF
03463 
03464 *----------------------------------------!
03465 * Table content
03466 * Note that by default RIGHT EDGE of bin is printed, as necessary for
03467 * cumulative distributions, this can be changed with SLAN option
03468 *----------------------------------------!
03469       CALL GLK_OptOut(idl(1),ioplog,iopsla,ioperb,iopsc1,iopsc2)
03470 *
03471 * table printout can be controlled by  GLK_SetTabRan(i1,i2,i3)
03472       k1=1
03473       k2=nchx
03474       k3=1
03475       IF( m_KeyTbr .EQ. 1 ) THEN
03476          k1 = MAX(k1,m_TabRan(1))
03477          k2 = MIN(k2,m_TabRan(2))
03478          k3 = MAX(k3,m_TabRan(3))
03479          m_KeyTbr = 0
03480       ENDIF
03481 *
03482       DO k=k1,k2,k3
03483          IF(chr2 .EQ. 'R') THEN
03484 * right
03485             xi= dxl + (dxu-dxl)*k/(1d0*nchx)
03486          ELSEIF(chr2 .EQ. 'L') THEN
03487 * left
03488             xi= dxl + (dxu-dxl)*(k-1d0)/(1d0*nchx)
03489          ELSE
03490 * midle
03491             xi= dxl + (dxu-dxl)*(k-0.5d0)/(1d0*nchx)
03492          ENDIF
03493          IF(ioperb.eq.2) THEN
03494           WRITE(m_ltx,'(A,'//fmt1//','//Cn(Npl)//'(A,'//fmt2//',A,A,'//fmt3//'),  A)')
03495      $                '$', xi, ('$ & $', bi(k,j), m_BS, 'pm', er(k,j), j=1,Npl), '$'
03496           WRITE(m_ltx,'(2A)') m_BS,m_BS
03497          ELSE
03498           WRITE(m_ltx,'(A,'//fmt1//','//Cn(Npl)//'(A,'//fmt2//'),  A)')
03499      $                '$', xi, ('$ & $', bi(k,j), j=1,Npl), '$'
03500           WRITE(m_ltx,'(2A)') m_BS,m_BS
03501          ENDIF
03502       ENDDO
03503 *------------------------------!
03504 * Ending
03505 *------------------------------!
03506       WRITE(m_ltx,'(4A)') m_BS,'hline',m_BS,'hline'
03507       WRITE(m_ltx,'(2A)') m_BS,'end{tabular}'
03508       IF(ABS(m_lint) .EQ. 2 ) THEN
03509          WRITE(m_ltx,'(A)') '% ========================================'
03510       ELSE
03511          WRITE(m_ltx,'(2A)') m_BS,'end{table}'
03512       ENDIF
03513       WRITE(m_ltx,'(A)') '% ============= end   table =============='
03514       WRITE(m_ltx,'(A)') '% ========================================'
03515       IF(chr3 .EQ. 'E') THEN
03516          WRITE(m_ltx,'(2A)') m_BS,'newpage'
03517       ELSE
03518          WRITE(m_ltx,'(A)') '% ========================================'
03519       ENDIF
03520       RETURN
03521  900  CALL GLK_Retu1(' ++++ GLK_PlTable2: Nonexistig histo,id= ',ID)
03522       RETURN
03523  901  CALL GLK_Retu1(' ++++ GLK_PlTable2: To many columns Nplt= ',Nplt)
03524       END
03525 
03526 
03527       SUBROUTINE GLK_WtMon(mode,id,par1,par2,par3)
03528 *     ********************************************
03529 * !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
03530 * !!!!  It is now replaces by GKL_M package, see below  !!!
03531 * !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
03532 *
03533 * Utility program for monitoring M.C. rejection weights.
03534 * ---------------------------------------------------------
03535 * It is backward compatible with WMONIT except:
03536 *  (1) for id=-1 one  should call as follows:
03537 *      GLK_WtMon(-1,id,0d0,1d0,1d0) or skip initialisation completely!
03538 *  (2) maximum absolute weight is looked for,
03539 *  (3) GLK_Print(-id) prints weight distribution, net profit!
03540 *  (4) no restriction id<100 any more!
03541 * ---------------------------------------------------------
03542 * wt is weight, wtmax is maximum weight and rn is random number.
03543 * IF(mode .EQ. -1) then
03544 *          initalization if entry id,
03545 *        - wtmax is maximum weight used for couting overweighted
03546 *          other arguments are ignored
03547 * ELSEIF(mode .EQ. 0) then
03548 *          summing up weights etc. for a given event for entry id
03549 *        - wt is current weight.
03550 *        - wtmax is maximum weight used for couting overweighted
03551 *          events with wt>wtmax.
03552 *        - rn is random number used in rejection, it is used to
03553 *          count no. of accepted (rn < wt/wtmax) and rejected
03554 *          (wt > wt/wtmax) events,
03555 *          if ro rejection then put rn=0d0.
03556 * ELSEIF(mode .EQ. 1) THEN
03557 *          in this mode wmonit repports on accumulated statistics
03558 *        - averwt= average weight wt counting all event
03559 *        - errela= relative error of averwt
03560 *        - nevtot= total number of accounted events
03561 *        - nevacc= no. of accepted events (rn < wt/wtmax)
03562 *        - nevneg= no. of events with negative weight (wt < 0)
03563 *        - nevzer= no. of events with zero weight (wt = 0d0)
03564 *        - nevove= no. of overweghted events (wt > wtmax)
03565 *          and if you do not want to use cmonit then the value
03566 *          the value of averwt is assigned to wt,
03567 *          the value of errela is assigned to wtmax and
03568 *          the value of wtmax  is assigned to rn in this mode.
03569 * ELSEIF(mode .EQ. 2) THEN
03570 *          all information defined for entry id defined above
03571 *          for mode=2 is just printed of unit nout
03572 * ENDIF
03573 * note that output repport (mode=1,2) is done dynamically just for a
03574 * given entry id only and it may be repeated many times for one id and
03575 * for various id's as well.
03576 *     ************************
03577       IMPLICIT NONE
03578       INCLUDE 'GLK.h'
03579       INTEGER           mode,id
03580       DOUBLE PRECISION  par1,par2,par3
03581 * locals
03582       INTEGER           idg,nevneg,nevzer,nevtot,nevove,nevacc,nbin,lact,ist3,ntot,ist,ist2
03583       DOUBLE PRECISION  xl,xu,errela,sswt,averwt,wwmax,swt,wt,wtmax,rn
03584 *---------------------------------------------------------------------------
03585       idg = -id
03586       IF(id .LE. 0) THEN
03587            CALL GLK_Stop1(' =====> GLK_WtMon: wrong id= ',id)
03588       ENDIF
03589       IF(mode .EQ. -1) THEN
03590 *     *******************
03591            nbin = nint(dabs(par3))
03592            IF(nbin .GT. 100) nbin =100
03593            IF(nbin .EQ. 0)   nbin =1
03594            xl   =  par1
03595            xu   =  par2
03596            IF(xu .LE. xl) THEN
03597              xl = 0d0
03598              xu = 1d0
03599            ENDIF
03600            CALL GLK_hadres(idg,lact)
03601            IF(lact .EQ. 0) THEN
03602               CALL GLK_Book1(idg,' GLK_WtMon $',nbin,xl,xu)
03603            ELSE
03604               WRITE(m_out,*) ' WARNING GLK_WtMon: exists, id= ',id
03605               WRITE(    6,*) ' WARNING GLK_WtMon: exists, id= ',id
03606            ENDIF
03607       ELSEIF(mode .EQ. 0) THEN
03608 *     **********************
03609          CALL GLK_hadres(idg,lact)
03610            IF(lact .EQ. 0) THEN
03611               WRITE(m_out,*) ' *****> GLK_WtMon: uninitialized, id= ',id
03612               WRITE(    6,*) ' *****> GLK_WtMon: uninitialized, id= ',id
03613               CALL GLK_Book1(idg,' GLK_WtMon $',1,0d0,1d0)
03614               CALL GLK_hadres(idg,lact)
03615            ENDIF
03616            wt   =par1
03617            wtmax=par2
03618            rn   =par3
03619 *     standard entries
03620            CALL GLK_Fil1(idg,wt,1d0)  !!!! <-- principal filling!!!!
03621 *     additional goodies
03622            ist  = m_index(lact,2)
03623            ist2 = ist+7
03624            ist3 = ist+11
03625 *    maximum weight -- maximum by absolute value but keeping sign
03626            m_b(ist3+13)    = max( dabs(m_b(ist3+13)) ,dabs(wt))
03627            IF(wt .NE. 0d0) m_b(ist3+13)=m_b(ist3+13) *wt/dabs(wt)
03628 *    nevzer,nevove,nevacc
03629            IF(wt .EQ. 0d0)        m_b(ist3+10) =m_b(ist3+10) +1d0
03630            IF(wt .GT. wtmax)      m_b(ist3+11) =m_b(ist3+11) +1d0
03631            IF(rn*wtmax .LE. wt)   m_b(ist3+12) =m_b(ist3+12) +1d0
03632       ELSEIF(mode .GE. 1 .OR. mode .LE. 10) THEN
03633 *     *************************************
03634          CALL GLK_hadres(idg,lact)
03635            IF(lact .EQ. 0) THEN
03636               CALL GLK_Stop1(' lack of initialization, id=',id)
03637            ENDIF
03638            ist    = m_index(lact,2)
03639            ist2   = ist+7
03640            ist3   = ist+11
03641            ntot   = nint(m_b(ist3 +7))
03642            swt    =      m_b(ist3 +8)
03643            sswt   =      m_b(ist3 +9)
03644            IF(ntot.LE.0 .OR. swt.EQ.0d0 )  THEN
03645               averwt=0d0
03646               errela=0d0
03647            ELSE
03648               averwt=swt/float(ntot)
03649               errela=sqrt(abs(sswt/swt**2-1d0/float(ntot)))
03650            ENDIF
03651            nevneg = m_b(ist3  +1) !!! it us underflow, xlow=0 assumed!!!
03652            nevzer = m_b(ist3 +10)
03653            nevove = m_b(ist3 +11)
03654            nevacc = m_b(ist3 +12)
03655            wwmax  = m_b(ist3 +13)
03656            nevtot = ntot
03657 * Output through parameters
03658            par1   = averwt
03659            par2   = errela
03660            par3   = nevtot
03661            IF(mode .EQ. 2) THEN
03662               par1   = nevacc
03663               par2   = nevneg
03664               par3   = nevove
03665            ELSEIF(mode .EQ. 3) THEN
03666               par1   = nevneg
03667               par2   = wwmax
03668            ENDIF
03669 *  no printout for mode <10
03670 *  ************************
03671            IF(mode .LE. 9) RETURN
03672            WRITE(m_out,1003) id, averwt, errela, wwmax
03673            WRITE(m_out,1004) nevtot,nevacc,nevneg,nevove,nevzer
03674            IF(mode .LE. 10) RETURN
03675            CALL GLK_Print(idg)
03676       ELSE
03677 *     ****
03678            CALL GLK_Stop1('+++GLK_WtMon: wrong mode=',mode)
03679       ENDIF
03680 *     *****
03681  1003 FORMAT(
03682      $  ' ======================= GLK_WtMon ========================='
03683      $/,'   id           averwt         errela            wwmax'
03684      $/,    i5,           e17.7,         f15.9,           e17.7)
03685  1004 FORMAT(
03686      $  ' -----------------------------------------------------------'
03687      $/,'      nevtot      nevacc      nevneg      nevove      nevzer'
03688      $/,   5i12)
03689       END
03690 
03691       SUBROUTINE GLK_CumHis(IdGen,id1,id2)
03692 *     ************************************
03693 *///////////////////////////////////////////////////////////////////////////
03694 *//   Cumulates histogram content starting from UNDERFLOW                 //
03695 *//   and normalizes to the total x-section in NANOBARNS                  //
03696 *//   IdGen is ID of special histogram written by M.C. generator itself   //
03697 *//   id2. NE. id1 required!!!                                            //
03698 *///////////////////////////////////////////////////////////////////////////
03699 *     ***********************************
03700       IMPLICIT NONE
03701       INTEGER  IdGen,id1,id2
03702 *----------------------------------------------------------------------
03703       INCLUDE 'GLK.h'
03704       SAVE
03705 *----------------------------------------------------------------------
03706       CHARACTER*80 TITLE
03707       DOUBLE PRECISION    X(m_MaxNb),ER(m_MaxNb)
03708       LOGICAL GLK_EXIST
03709       DOUBLE PRECISION    swt,sswt,xsec,errel,tmin,tmax
03710       DOUBLE PRECISION    xscrnb,ERela,WtSup
03711       INTEGER  i,nbt,nevt
03712       DOUBLE PRECISION    GLK_hi,GLK_hie
03713 *----------------------------------------------------------------------
03714       IF (GLK_Exist(id2)) GOTO 900
03715 *
03716       CALL GLK_MgetNtot(IdGen,nevt)
03717       CALL GLK_MgetAve( IdGen,xscrnb,ERela,WtSup)
03718 *
03719       IF(nevt .EQ. 0) GOTO 901
03720       CALL GLK_hinbo1(id1,title,nbt,tmin,tmax)
03721       swt  = GLK_hi( id1,0)     ! UNDERFLOW
03722       sswt = GLK_hie(id1,0)**2  ! UNDERFLOW
03723       DO i=1,nbt
03724          swt   = swt + GLK_hi( id1,i)
03725          sswt  = sswt+ GLK_hie(id1,i)**2
03726 * note NEVT in error calc. is for the entire sample related
03727 * to the crude x-section XCRU including !!! zero weight events !!!!
03728          xsec  = 0d0
03729          errel = 0d0
03730          IF(swt .NE. 0d0 .AND. nevt .NE. 0) THEN
03731             xsec  = swt*(xscrnb/nevt)
03732             errel = SQRT(ABS(sswt/swt**2-1d0/FLOAT(nevt)))
03733          ENDIF
03734          x(i)  = xsec
03735          er(i) = xsec*errel
03736       ENDDO
03737 *! store result in id2
03738       CALL GLK_Book1(id2,title,nbt,tmin,tmax)
03739       CALL GLK_Pak(  id2,x)
03740       CALL GLK_Pake( id2,er)
03741       CALL GLK_idopt(id2,'ERRO')
03742       RETURN
03743  900  WRITE(6,*) '+++++ CUMHIS: ID2 exixsts!!',ID2
03744       RETURN
03745  901  WRITE(6,*) '+++++ CUMHIS: EMPTY HISTO ID=',ID1
03746       END
03747 
03748 
03749 
03750 
03751       SUBROUTINE GLK_RenHst(chak,IdGen,id1,id2)
03752 *     *****************************************
03753 *///////////////////////////////////////////////////////////////////////////
03754 *//   IdGen is ID of special histogram written by M.C. generator itself   //
03755 *//   This routine RE-NORMALIZES to  NANOBARNS or to UNITY                //
03756 *//   CHAK = 'NB  '    normal case [nb]                                   //
03757 *//   CHAK = 'NB10'    log10 x-scale assumed [nb]                         //
03758 *//   CHAK = 'UNIT'    normalization to unity                             //
03759 *//   id2 .NE. id1 required !!!                                           //
03760 *///////////////////////////////////////////////////////////////////////////
03761 *     ***********************************
03762       IMPLICIT NONE
03763       CHARACTER*4 CHAK
03764       INTEGER     IdGen,id1,id2
03765 *----------------------------------------------------------------------
03766       INCLUDE 'GLK.h'
03767       SAVE
03768       CHARACTER*80 TITLE
03769       DOUBLE PRECISION        xscrnb,ERela,WtSup,tmin,tmax
03770       DOUBLE PRECISION        swt,fln10,fact
03771       INTEGER      i,nbt,nevt
03772       DOUBLE PRECISION    GLK_hi,GLK_hie
03773 *----------------------------------------------------------------------
03774       IF( id2 .eq. id1) GOTO 900
03775 
03776       CALL GLK_MgetNtot(IdGen,nevt)
03777       CALL GLK_MgetAve( IdGen,xscrnb,ERela,WtSup)
03778 *
03779       CALL GLK_hinbo1(id1,title,nbt,tmin,tmax)
03780       IF(     chak .EQ. 'NB  ') THEN
03781          fact = nbt*xscrnb/(nevt*(tmax-tmin))
03782          CALL GLK_Operat(id1,'+',id1,id2, fact, 0d0)
03783       ELSEIF( chak .EQ. 'NB10') THEN
03784          fln10 = log(10.)
03785          fact = nbt*xscrnb/(nevt*(tmax-tmin)*fln10)
03786          CALL GLK_Operat(id1,'+',id1,id2, fact, 0d0)
03787       ELSEIF( chak .EQ. 'UNIT') THEN
03788          swt  = GLK_hi(id1,0)
03789          DO i=1,nbt+1
03790             swt   = swt + GLK_hi(id1,i)
03791          ENDDO
03792          fact = nbt/((tmax-tmin))/swt
03793          CALL GLK_Operat(id1,'+',id1,id2, fact, 0d0)
03794       ELSEIF( chak .EQ. 'UN10') THEN
03795          swt  = GLK_hi(id1,0)
03796          DO i=1,nbt+1
03797             swt   = swt + GLK_hi(id1,i)
03798          ENDDO
03799          fact = nbt/((tmax-tmin)*log(10.))/swt
03800          CALL GLK_Operat(id1,'+',id1,id2, fact, 0d0)
03801       ELSEIF( chak .EQ. '    ') THEN
03802          CALL GLK_Operat(id1,'+',id1,id2, 1d0, 0d0)
03803       ELSE
03804          WRITE(6,*) '+++++ RENHST: wrong chak=',chak
03805       ENDIF
03806 *
03807       RETURN
03808  900  WRITE(6,*) '+++++ RENHST: ID1=ID2=',ID1
03809       END
03810 
03811 
03812 *///////////////////////////////////////////////////////////////////////
03813 *//                New Weight Motoring ToolBox
03814 *//               (replacement for WTmonit etc.)
03815 *//
03816 *//  The tool to monitor very precisely the average weigh
03817 *//  and other features of the weight distribution.
03818 *//  Note that in principle we are vitaly interested in three parts
03819 *//  of the weight distribution:
03820 *//     Underflow (-infty,0)
03821 *//     Regular   (0, WTmax)
03822 *//     Overflow  (WTmax,+infty)
03823 *//  with special emphasis on events with exactly zero weight WT=0d0.
03824 *//  Nevertheless, we split (0, WTmax) range into several bins
03825 *//  in order to be able to visualise the weight distribution.
03826 *//  (Using stardard tools for histogram)
03827 *//
03828 *//
03829 *///////////////////////////////////////////////////////////////////////
03830       SUBROUTINE GLK_Mbook(idm,title,nnchx,WTmax)
03831 *     ******************************************
03832 *///////////////////////////////////////////////////////////////////////
03833 *//
03834 *//   Booking one entry. Note it is not an ordinary histogram!!!
03835 *//   It works just like GLK_Book1 except that it
03836 *//   has internaly negative id and x_minimum is always zero.
03837 *//
03838 *///////////////////////////////////////////////////////////////////////
03839       IMPLICIT NONE
03840       INCLUDE 'GLK.h'
03841       SAVE
03842       INTEGER idm
03843       CHARACTER*80 title
03844       DOUBLE PRECISION   WTmax
03845 *
03846       LOGICAL GLK_Exist
03847       INTEGER j,id,nnchx,nchx,lact,lengt2,ist,ist2,ist3
03848       INTEGER iopsc1, iopsc2, ioperb, ioplog, iopsla
03849       INTEGER iflag1, iflag2
03850       INTEGER ityphi
03851       DOUBLE PRECISION   xl,xu,ddx
03852 *-------------------------------------------------
03853       CALL GLK_Initialize
03854       id = -idm
03855       IF(GLK_Exist(id)) goto 900
03856       ist=m_length
03857       CALL GLK_hadres(0,lact)
03858 * Check if there is a free entry in the m_index
03859       IF(lact .EQ. 0) CALL GLK_Stop1('GLK_Mbook: no space left,id= ',id)
03860       m_index(lact,1)=id
03861       m_index(lact,2)=m_length
03862       m_index(lact,3)=0
03863 * ---------- limits 
03864       CALL GLK_Copch(title,m_titlc(lact))
03865       nchx =nnchx
03866       IF(nchx .GT. m_MaxNb)
03867      $     CALL GLK_Stop1(' GLK_Mbook: Too many bins ,id= ',id)
03868       xl   = 0d0
03869       xu   = WTmax
03870 * ---------- title and bin content ----------
03871       lengt2 = m_length +2*nchx +m_buf1+1
03872       IF(lengt2 .GE. m_LenmB)
03873      $  CALL GLK_Stop1('GLK_Mbook:too litle storage, m_LenmB= ',m_LenmB)
03874 *
03875       DO j=m_length+1,lengt2+1
03876          m_b(j) = 0d0
03877       ENDDO
03878       m_length=lengt2
03879 *... default flags
03880       ioplog   = 1
03881       iopsla   = 1
03882       ioperb   = 1
03883       iopsc1   = 1
03884       iopsc2   = 1
03885       iflag1   =
03886      $ ioplog+10*iopsla+100*ioperb+1000*iopsc1+10000*iopsc2
03887       ityphi   = 3  !!!! <-- new type of histo !!!!
03888       iflag2   = ityphi
03889 * examples of decoding flags
03890 *      id       = nint(m_b(ist+2)-9d0-9d12)/10
03891 *      iflag1   = nint(m_b(ist+3)-9d0-9d12)/10
03892 *      ioplog = mod(iflag1,10)
03893 *      iopsla = mod(iflag1,100)/10
03894 *      ioperb = mod(iflag1,1000)/100
03895 *      iopsc1 = mod(iflag1,10000)/1000
03896 *      iopsc2 = mod(iflag1,100000)/10000
03897 *      iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
03898 *      ityphi = mod(iflag2,10)
03899 *--------- buffer -----------------
03900 * header
03901       m_b(ist +1)  = 9999999999999d0
03902       m_b(ist +2)  = 9d12 +     id*10 +9d0
03903       m_b(ist +3)  = 9d12 + iflag1*10 +9d0
03904       m_b(ist +4)  = 9d12 + iflag2*10 +9d0
03905 * dummy vertical scale
03906       m_b(ist +5)  =  -100d0
03907       m_b(ist +6)  =   100d0
03908 * pointer used to speed up search of histogram address
03909       m_b(ist +7)  =   0d0
03910 * information on binning
03911       ist2         = ist+7
03912       m_b(ist2 +1) = nchx
03913       m_b(ist2 +2) = xl
03914       m_b(ist2 +3) = xu
03915       ddx = xu-xl
03916       IF(ddx .EQ. 0d0)
03917      $     CALL GLK_Stop1(' GLK_Mbook:    xl=xu,              id= ',id)
03918       m_b(ist2 +4) = float(nchx)/ddx
03919 *
03920 * under/over-flow etc.
03921       ist3       = ist+11
03922       DO j=1,13
03923          m_b(ist3 +j)=0d0
03924       ENDDO
03925       RETURN
03926 *----------------
03927  900  CALL GLK_Retu1(' WARNING GLK_Mbook: already exists id= ', id)
03928       END
03929 
03930 
03931       SUBROUTINE GLK_Mfill(idm,Wtm,rn)
03932 *     ********************************
03933 *///////////////////////////////////////////////////////////////////////
03934 *//
03935 *//    filling of M-subpackage entry
03936 *//    simillar as fil1 for 1-dim histo but the storage for error
03937 *//    is now used to store sum for 'partial averages' <wt-xlowedge>
03938 *//
03939 *///////////////////////////////////////////////////////////////////////
03940       IMPLICIT NONE
03941       INTEGER idm
03942       DOUBLE PRECISION  Wtm,rn
03943       INCLUDE 'GLK.h'
03944       SAVE
03945       INTEGER id
03946       INTEGER lact, ist, ist2, ist3, iflag2, nchx, ityphi
03947       INTEGER iposx1,ipose1, kposx1, kpose1, kx
03948       DOUBLE PRECISION   Wt, deltx, factx, xlowedge
03949       DOUBLE PRECISION   xu, xl, x1, wtmax
03950 *---------------------------------
03951       id = -idm
03952       Wt = Wtm
03953       CALL GLK_hadres(id,lact)
03954 * exit for non-existig histo
03955       IF(lact .EQ. 0)  
03956      $     CALL GLK_Stop1('+++GLK_Mfill: nonexisting id= ',id)
03957 
03958       ist  = m_index(lact,2)
03959       ist2 = ist+7
03960       ist3 = ist+11
03961 * one-dim. histo only
03962       iflag2   = nint(m_b(ist+4)-9d0-9d12)/10
03963       ityphi   = mod(iflag2,10)
03964       IF(ityphi .NE. 3) CALL GLK_Stop1('+++GLK_Mfill: wrong id=  ',id)
03965       x1 =  Wt
03966       m_index(lact,3)=m_index(lact,3)+1
03967 * for standard average of x=Wt and its error
03968       m_b(ist3 +7)  =m_b(ist3 +7)  +1
03969       m_b(ist3 +8)  =m_b(ist3 +8)  +x1
03970       m_b(ist3 +9)  =m_b(ist3 +9)  +x1*x1
03971 * filling coordinates
03972       nchx  = m_b(ist2 +1)
03973       xl    = m_b(ist2 +2)      !!<--- It was set to zero in book!!!
03974       xu    = m_b(ist2 +3)
03975       WtMax = xu
03976       factx = m_b(ist2 +4)      ! (fact=nchx/(xu-xl)
03977       deltx = 1d0/factx
03978       IF(x1 .LT. xl) THEN
03979 * (U)nderflow
03980          iposx1 = ist3 +1
03981          ipose1 = ist3 +4
03982          m_b(iposx1) = m_b(iposx1)  +1d0
03983          m_b(ipose1) = m_b(ipose1)  +Wt
03984       ELSEIF(x1 .GT. xu) THEN
03985 * (O)verflow
03986          iposx1 = ist3 +3
03987          ipose1 = ist3 +6
03988          kposx1 = 0
03989          m_b(iposx1) = m_b(iposx1)  +1d0
03990          m_b(ipose1) = m_b(ipose1)  +(Wt- WtMax)
03991       ELSE
03992 * All of (R)egular range (0,WtMax) in one bin
03993          iposx1 = ist3 +2
03994          ipose1 = ist3 +5
03995          m_b(iposx1) = m_b(iposx1)  +1d0
03996          m_b(ipose1) = m_b(ipose1)  +Wt
03997 * (R)egular bin, the ACTUAL one
03998          kx = (x1-xl)*factx+1d0
03999          kx = MIN( MAX(kx,1) ,nchx)
04000          kposx1 = ist +m_buf1+kx
04001          kpose1 = ist +m_buf1+nchx+kx
04002          xlowedge = deltx*(kx-1)
04003          m_b(kposx1) = m_b(kposx1)  +1d0
04004          m_b(kpose1) = m_b(kpose1)  +(Wt-xlowedge)
04005       ENDIF
04006 *--------------------------------
04007 * Additional goodies:
04008 * maximum weight -- maximum by absolute value but keeping sign
04009       m_b(ist3+13) = MAX( DABS(m_b(ist3+13)) ,DABS(wt))
04010       IF(wt .NE. 0d0) m_b(ist3+13)=m_b(ist3+13) *wt/dabs(wt)
04011 * nevzer,nevove,nevacc
04012       IF(wt .EQ. 0d0)        m_b(ist3+10) =m_b(ist3+10) +1d0
04013       IF(wt .GT. wtmax)      m_b(ist3+11) =m_b(ist3+11) +1d0
04014       IF(rn*wtmax .LE. wt)   m_b(ist3+12) =m_b(ist3+12) +1d0
04015 *---
04016       END   !GLK_Mfill
04017 
04018 
04019       SUBROUTINE GLK_MgetAll(idm,
04020      $     AveWt,ERela, WtSup, AvUnd, AvOve,
04021      $     Ntot,Nacc,Nneg,Nove,Nzer)
04022 *     ***************************************************************
04023 *///////////////////////////////////////////////////////////////////////
04024 *//
04025 *//   Get all statistics out
04026 *//
04027 *///////////////////////////////////////////////////////////////////////
04028       IMPLICIT NONE
04029       INTEGER idm
04030       DOUBLE PRECISION   AveWt,ERela, WtSup, AvUnd, AvOve
04031       INTEGER Ntot,Nacc,Nneg,Nove,Nzer
04032       INCLUDE 'GLK.h'
04033       SAVE
04034       INTEGER id,ist,ist2,ist3,lact
04035       DOUBLE PRECISION   swt,sswt
04036 *--------------------
04037       id= -idm
04038       CALL GLK_hadres(id,lact)
04039       IF(lact .EQ. 0)
04040      $     CALL GLK_Stop1('GLK_MgetAll:lack of initialization, id=',id)
04041       ist    = m_index(lact,2)
04042       ist2   = ist+7
04043       ist3   = ist+11
04044       Ntot   = nint(m_b(ist3 +7))
04045       swt    =      m_b(ist3 +8)
04046       sswt   =      m_b(ist3 +9)
04047       IF(Ntot.LE.0 .OR. swt.EQ.0d0 )  THEN
04048          AveWt=0d0
04049          ERela=0d0
04050       ELSE
04051          AveWt=swt/DFLOAT(Ntot)
04052          ERela=sqrt(abs(sswt/swt**2-1d0/float(Ntot)))
04053       ENDIF
04054       WtSup  = m_b(ist3 +13)
04055       AvUnd  = m_b(ist3 +4)/Ntot
04056       AvOve  = m_b(ist3 +6)/Ntot
04057       Nneg = m_b(ist3  +1)    ! NB. it is underflow
04058       Nzer = m_b(ist3 +10)
04059       Nove = m_b(ist3 +11)
04060       Nacc = m_b(ist3 +12)
04061 *-----------------------------
04062 *      WRITE(m_out,1003) idm, AveWt, ERela, WtSup
04063 *      WRITE(m_out,1004) Ntot,Nacc,Nneg,Nove,Nzer
04064 * 1003 FORMAT(
04065 *     $  ' ======================= GLK_Mget =========================='
04066 *     $/,'   id            AveWt          ERela            WtSup'
04067 *     $/,    i5,           e17.7,         f15.9,           e17.7)
04068 * 1004 FORMAT(
04069 *     $  ' -----------------------------------------------------------'
04070 *     $/,'        Ntot        Nacc        Nneg        Nove        Nzer'
04071 *     $/,   5i12)
04072 *------------------------------
04073       END
04074 
04075       SUBROUTINE GLK_MgetNtot(id,Ntot)
04076 *///////////////////////////////////////////////////////////////////////
04077 *//
04078 *//   Get Ntotal only
04079 *//
04080 *///////////////////////////////////////////////////////////////////////
04081       IMPLICIT NONE
04082       INCLUDE 'GLK.h'
04083       SAVE
04084       INTEGER idm,id
04085       DOUBLE PRECISION   AveWt, ERela, WtSup, AvUnd, AvOve
04086       INTEGER Ntot, Nacc, Nneg, Nove, Nzer
04087 *--------------------
04088       CALL GLK_MgetAll(id,
04089      $     AveWt,ERela, WtSup, AvUnd, AvOve,
04090      $     Ntot,Nacc,Nneg,Nove,Nzer)
04091       END
04092 
04093       SUBROUTINE GLK_MgetAve(id,AveWt,ERela,WtSup)
04094 *///////////////////////////////////////////////////////////////////////
04095 *//
04096 *//   Get averages only and highest weight
04097 *//
04098 *///////////////////////////////////////////////////////////////////////
04099       IMPLICIT NONE
04100       INCLUDE 'GLK.h'
04101       SAVE
04102       INTEGER idm,id
04103       DOUBLE PRECISION   AveWt, ERela, WtSup, AvUnd, AvOve
04104       INTEGER Ntot, Nacc, Nneg, Nove, Nzer
04105 *--------------------
04106       CALL GLK_MgetAll(id,
04107      $     AveWt,ERela, WtSup, AvUnd, AvOve,
04108      $     Ntot,Nacc,Nneg,Nove,Nzer)
04109       END
04110 
04111       SUBROUTINE GLK_Mprint(idm)
04112 *///////////////////////////////////////////////////////////////////////
04113 *//
04114 *//   Printout
04115 *//   Note that bin errors have now changed meaning
04116 *//
04117 *///////////////////////////////////////////////////////////////////////
04118       IMPLICIT NONE
04119       INCLUDE 'GLK.h'
04120       SAVE
04121       INTEGER idm,id
04122       id= -idm
04123       CALL GLK_Print(id)
04124       END
04125 
04126 *//////////////////////////////////////////////////////////////////////////////
04127 *//////////////////////////////////////////////////////////////////////////////
04128 *//////////////////////////////////////////////////////////////////////////////
04129 *//                                                                          //
04130 *//                     Setters and Getters                                  //
04131 *//                                                                          //
04132 *//////////////////////////////////////////////////////////////////////////////
04133 
04134       SUBROUTINE GLK_Clone1(id1,id2,title2)
04135 *//////////////////////////////////////////////////////////////////////////////
04136 *//   Clone 1-dimensional histo with onl bining and new title                //
04137 *//////////////////////////////////////////////////////////////////////////////
04138       CHARACTER*80 title1, title2, title3
04139       INTEGER      i,nb
04140       DOUBLE PRECISION        xmin,xmax
04141 
04142       CALL GLK_hinbo1(id1,title1,nb,xmin,xmax)
04143       CALL GLK_Copch(title2,title3)
04144       CALL GLK_Book1(id2,title3,nb,xmin,xmax)
04145 
04146       END
04147 
04148       SUBROUTINE GLK_Ymimax(id,wmin,wmax)
04149 *//////////////////////////////////////////////////////////////////////////////
04150 *//                                                                          //
04151 *//////////////////////////////////////////////////////////////////////////////
04152       IMPLICIT NONE
04153       INTEGER id
04154       DOUBLE PRECISION  wmin,wmax
04155 
04156       CALL GLK_Yminim(id,wmin)
04157       CALL GLK_Ymaxim(id,wmax)
04158       END
04159 
04160 
04161       SUBROUTINE GLK_PLset(ch,xx)
04162 *//////////////////////////////////////////////////////////////////////////////
04163 *//                                                                          //
04164 *//   Old style seter, sets type of the ploting mark                         //
04165 *//                                                                          //
04166 *//////////////////////////////////////////////////////////////////////////////
04167       IMPLICIT NONE
04168       CHARACTER*4 CH
04169       DOUBLE PRECISION  xx
04170       INCLUDE 'GLK.h'
04171       SAVE
04172 *----------------------------------
04173       IF(CH .EQ. 'DMOD') THEN
04174         m_tline = NINT(xx)
04175       ENDIF
04176       END
04177 
04178       SUBROUTINE GLK_SetNout(ilun)
04179 *//////////////////////////////////////////////////////////////////////////////
04180 *//                                                                          //
04181 *//   Sets output unit number                                                //
04182 *//                                                                          //
04183 *//////////////////////////////////////////////////////////////////////////////
04184       IMPLICIT NONE
04185       INCLUDE 'GLK.h'
04186       SAVE
04187       INTEGER ilun
04188 
04189       CALL GLK_Initialize
04190       m_out=ilun
04191       END
04192 
04193       SUBROUTINE GLK_GetYmin(id,ymin)
04194 *//////////////////////////////////////////////////////////////////////////////
04195 *//   Sets vertical scale                                                    //
04196 *//////////////////////////////////////////////////////////////////////////////
04197       IMPLICIT NONE
04198       INCLUDE 'GLK.h'
04199       INTEGER id
04200       DOUBLE PRECISION   ymin
04201       INTEGER lact,ist
04202 *
04203       CALL GLK_hadres(id,lact)
04204       IF(lact .EQ. 0) RETURN
04205       ist= m_index(lact,2)
04206       ymin = m_b(ist+5)
04207       END
04208 
04209       SUBROUTINE GLK_GetYmax(id,ymax)
04210 *//////////////////////////////////////////////////////////////////////////////
04211 *//   Sets vertical scale                                                    //
04212 *//////////////////////////////////////////////////////////////////////////////
04213       IMPLICIT NONE
04214       INCLUDE 'GLK.h'
04215       INTEGER id
04216       DOUBLE PRECISION   ymax
04217       INTEGER lact,ist
04218 *
04219       CALL GLK_hadres(id,lact)
04220       IF(lact .EQ. 0) RETURN
04221       ist= m_index(lact,2)
04222       ymax = m_b(ist+6)
04223       END
04224 
04225       SUBROUTINE GLK_SetYmin(id,ymin)
04226 *//////////////////////////////////////////////////////////////////////////////
04227 *//   Sets vertical scale                                                    //
04228 *//////////////////////////////////////////////////////////////////////////////
04229       IMPLICIT NONE
04230       INCLUDE 'GLK.h'
04231       INTEGER id
04232       DOUBLE PRECISION   ymin
04233       INTEGER lact,ist
04234 *
04235       CALL GLK_hadres(id,lact)
04236       IF(lact .EQ. 0) RETURN
04237       ist= m_index(lact,2)
04238       m_b(ist+5) = ymin
04239       CALL GLK_idopt(id,'YMIN')
04240       END
04241 
04242       SUBROUTINE GLK_SetYmax(id,ymax)
04243 *//////////////////////////////////////////////////////////////////////////////
04244 *//   Sets vertical scale                                                    //
04245 *//////////////////////////////////////////////////////////////////////////////
04246       IMPLICIT NONE
04247       INCLUDE 'GLK.h'
04248       INTEGER id
04249       DOUBLE PRECISION   ymax
04250       INTEGER lact,ist
04251 *
04252       CALL GLK_hadres(id,lact)
04253       IF(lact .EQ. 0) RETURN
04254       ist= m_index(lact,2)
04255       m_b(ist+6) = ymax
04256       CALL GLK_idopt(id,'YMAX')
04257       END
04258 
04259 
04260       SUBROUTINE GLK_GetYminYmax(id,ymin,ymax)
04261 *//////////////////////////////////////////////////////////////////////////////
04262 *//   Sets vertical scale                                                    //
04263 *//////////////////////////////////////////////////////////////////////////////
04264       IMPLICIT NONE
04265       INCLUDE 'GLK.h'
04266       INTEGER id
04267       DOUBLE PRECISION   ymin,ymax
04268 *
04269       CALL GLK_GetYmin(id,ymin)
04270       CALL GLK_GetYmax(id,ymax)
04271       END
04272 
04273       SUBROUTINE GLK_SetYminYmax(id,ymin,ymax)
04274 *//////////////////////////////////////////////////////////////////////////////
04275 *//   Sets vertical scale                                                    //
04276 *//////////////////////////////////////////////////////////////////////////////
04277       IMPLICIT NONE
04278       INCLUDE 'GLK.h'
04279       INTEGER id
04280       DOUBLE PRECISION   ymin,ymax
04281 *
04282       CALL GLK_SetYmin(id,ymin)
04283       CALL GLK_SetYmax(id,ymax)
04284       END
04285 
04286       SUBROUTINE GLK_CopyYmin(id1,id2)
04287 *//////////////////////////////////////////////////////////////////////////////
04288 *//   Sets vertical scale                                                    //
04289 *//////////////////////////////////////////////////////////////////////////////
04290       IMPLICIT NONE
04291       INCLUDE 'GLK.h'
04292       INTEGER id1,id2
04293       DOUBLE PRECISION   ymin
04294 *
04295       CALL GLK_GetYmin(id1,ymin)
04296       CALL GLK_SetYmin(id2,ymin)
04297       END
04298 
04299       SUBROUTINE GLK_CopyYmax(id1,id2)
04300 *//////////////////////////////////////////////////////////////////////////////
04301 *//   Sets vertical scale                                                    //
04302 *//////////////////////////////////////////////////////////////////////////////
04303       IMPLICIT NONE
04304       INCLUDE 'GLK.h'
04305       INTEGER id1,id2
04306       DOUBLE PRECISION   ymax
04307 *
04308       CALL GLK_GetYmax(id1,ymax)
04309       CALL GLK_SetYmax(id2,ymax)
04310       END
04311 
04312       SUBROUTINE GLK_SetColor(Color)
04313 *//////////////////////////////////////////////////////////////////////////////
04314 *//                                                                          //
04315 *//   Sets output unit number                                                //
04316 *//                                                                          //
04317 *//////////////////////////////////////////////////////////////////////////////
04318       IMPLICIT NONE
04319       INCLUDE 'GLK.h'
04320       CHARACTER*(*) Color
04321 *
04322       CALL GLK_Copch(Color,m_Color)
04323 *
04324       m_KeyCol = 1              !flag set up
04325       END
04326 
04327       SUBROUTINE GLK_SetTabRan(i1,i2,i3)
04328 *//////////////////////////////////////////////////////////////////////////////
04329 *//                                                                          //
04330 *//   Sets table range for taple printout in GKL_PlTable2                    //
04331 *//   i1,i2,i3 are lower limit, upper limit and increment                    //
04332 *//                                                                          //
04333 *//////////////////////////////////////////////////////////////////////////////
04334       IMPLICIT NONE
04335       INCLUDE 'GLK.h'
04336       INTEGER i1,i2,i3
04337 *
04338       m_KeyTbr = 1              !flag set up
04339       m_TabRan(1) = i1
04340       m_TabRan(2) = i2
04341       m_TabRan(3) = i3
04342       END
04343 
04344       SUBROUTINE GLK_GetNb(id,Nb)
04345 *//////////////////////////////////////////////////////////////////////////////
04346 *//                                                                          //
04347 *//////////////////////////////////////////////////////////////////////////////
04348       IMPLICIT NONE
04349       INCLUDE 'GLK.h'
04350       INTEGER  id,Nb
04351 * local
04352       CHARACTER*80 title
04353       INTEGER lact,ist,ist2
04354 *-------------
04355       CALL GLK_hadres(id,lact)
04356       IF(lact .EQ. 0) THEN
04357          CALL GLK_Stop1('+++STOP in GLK_GetNb: wrong id=',id)
04358       ENDIF
04359       ist    = m_index(lact,2)
04360       ist2   = ist+7
04361       Nb     = m_b(ist2 +1)
04362       END
04363 
04364       SUBROUTINE GLK_GetBin(id,ib,Bin)
04365 *//////////////////////////////////////////////////////////////////////////////
04366 *//                                                                          //
04367 *//  getting out bin content                                                 //
04368 *//                                                                          //
04369 *//////////////////////////////////////////////////////////////////////////////
04370       IMPLICIT NONE
04371       INCLUDE 'GLK.h'
04372       INTEGER    id,ib
04373       DOUBLE PRECISION      Bin,GLK_hi
04374 
04375       Bin =  GLK_hi(id,ib)
04376       END
04377 
04378 
04379 *//////////////////////////////////////////////////////////////////////////////
04380 *//                                                                          //
04381 *//                      End of CLASS  GLK                                   //
04382 *//////////////////////////////////////////////////////////////////////////////
04383 
Generated on Sun Oct 20 20:24:08 2013 for C++InterfacetoTauola by  doxygen 1.6.3