MUON/algo.F

   1       SUBROUTINE trig(y,x,iflag)
   2 C
   3 C *** DIGITISATION FOR THE MUON RAW DATA AFTER EACH EVENT ***
   4 C *** NVE 24-SEP-1990 CERN GENEVA ***
   5 C
   6 C CALLED BY : GUDIGI
   7 C ORIGIN    : NICK VAN EIJNDHOVEN
   8 C
   9 c
  10 c version 2 (open geom + L.U.T. + option DS on)
  11 c Input  : Hits on the muon trigger chambers (RPCs) from GALICE
  12 c Output : muon trigger decision (Unlike Sign) L0 low and high Pt
  13 c
  14 c 1- TRIGMAP --> DESCRIBE MUON trigger GEOMETRY
  15 c 2- REMPL   --> FILL bit pattern of MUON trigger CIRCUITS
  16 c 3- ALGO    --> INDICATE MUON trigger DECISION
  17 c
  18 c
  19
  20       real x(4,1000),y(4,1000)     ! 1000=nhitmax
  21       common/debug/idebug
  22
  23       data nhitmax/1000/           ! max possible number of hits
  24       data idebug/0/               ! for debuging
  25
  26 c      call hropen(99,'hist','hist_paw.hbook',' ',1024,ISTAT)
  27 c      open(UNIT=99,file='hist_paw.hbook',recl=1024,form='unformatted',
  28 c     +   access='direct',status='unknown')
  29 c      call hlimit(500000)
  30 c      call hbook2(1,'X-Y hits plane 1(1600)',100,-300.,300.,
  31 c     +                                       100,-300.,300.,0.)
  32 c      call hbook2(2,'X-Y hits plane 2(1615)',100,-300.,300.,
  33 c     +                                       100,-300.,300.,0.)
  34 c      call hbook2(3,'X-Y hits plane 3(1700)',100,-300.,300.,
  35 c     +                                       100,-300.,300.,0.)
  36 c      call hbook2(4,'X-Y hits plane 4(1715)',100,-300.,300.,
  37 c     +                                       100,-300.,300.,0.)
  38 cc
  39 c      call hbook1(10,'NUM of non-empty circuits',346,0.,346.,0.)
  40 c      call hbook1(11,'DISTR of non-empty circuits',
  41 c     +                346,-173.,173.,0.)
  42 c      call hbook1(12,'NUM of non-empty X circuits',346,0.,346.,0.)
  43 c      call hbook1(13,'NUM of non-empty Y circuits',346,0.,346.,0.)
  44 cc
  45 c
  46 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  47 c Describe trigger geometry
  48       call TRIGMAPP
  49 c
  50 c Hits on trigger chambers (RPCs) from GALICE are passed through
  51 c      x(iplan,ihit) and y(iplan,ihit)
  52 c with iplan=1-4 corresp. to Z=1600,1615,1700,1715 cm respectively
  53 c      ihit =number associated to a hit (ihit<1000, see nhitmax)
  54 c
  55 c
  56 c init
  57 c        do iplan=1,4
  58 c         do ihit=1,nhitmax
  59 c           x(iplan,ihit)=0
  60 c           y(iplan,ihit)=0
  61 c         enddo
  62 c        enddo
  63 c
  64 c for debugging proposes : ex of test sequences
  65 c        call seqtest
  66
  67 c for debugging proposes : test full evt pbpb simulated with geant
  68 c        call FULLEVT
  69
  70 c
  71 c Associate Hits to circuits (bit pattern)
  72             call REMPL(x,y,nhitmax)
  73 c Look for possible triggers (L0 , Low and High Pt cut at present)
  74             call ALGO(itrigR,itrigL0,itrigH0)
  75
  76 c  TRIGGER ?
  77
  78        print *,' '
  79        print *,'DIMUON TRIGGER - DIMUON TRIGGER'
  80        print *,'ROAD +/-8   ',itrigR
  81        print *,'Low  pt L0  ',itrigL0
  82        print *,'High pt L0  ',itrigH0
  83
  84 c      call hrfile(99,'hist','n')
  85 c      call hrout(0,icycle,' ')
  86 c      call hrend('hist')
  87 c      stop
  88        iflag=100*itrigH0+10*itrigL0+itrigR
  89       end
  90 c
  91 c
  92 c
  93       subroutine TRIGMAPP
  94 c
  95 c to be called at a begining of a run for dimuon trigger GEOM MAPing
  96 c
  97 c
  98 c      integer NUM(173),CODEX(173),CODEY(173)           !173=Ncirc
  99 c      real xcmax(4,173),xcmin(4,173),
 100 c     +     ycmax(4,173),ycmin(4,173),
 101 c     +     xycmax(4,173),xycmin(4,173)
 102 c      real xwi_c(173),xwi_m(173),xwi_p(173),ywi(173)
 103 c      integer nstrip_c(173),nstrip_p(173),nstrip_m(173)
 104 c
 105 c      common/TRIGMAP/xcmax,xcmin,ycmax,ycmin,xycmax,xycmin,
 106 c     +               CODEX,CODEY,
 107 c     +               xwi_c,xwi_m,xwi_p,ywi,
 108 c     +               nstrip_c,nstrip_m,nstrip_p,
 109 c     +               NUM,Z1,Z2,Z3,Z4,Ncirc
 110
 111       common/TRIGMAP/xcmax(4,173),xcmin(4,173),ycmax(4,173),
 112      +               ycmin(4,173),xycmax(4,173),xycmin(4,173),
 113      +               CODEX(173),CODEY(173),
 114      +               xwi_c(173),xwi_m(173),xwi_p(173),ywi(173),
 115      +               nstrip_c(173),nstrip_m(173),nstrip_p(173),
 116      +               NUM(173),Z1,Z2,Z3,Z4,Ncirc
 117       integer nstrip_c,nstrip_p,nstrip_m
 118       integer NUM,CODEX,CODEY,Ncirc
 119
 120       common/debug/idebug
 121
 122       integer ixwi_c,ixwi_p,ixwi_m
 123       integer iywi,iywi_c,iywi_p,iywi_m,ixyco_m,ixyco_p
 124
 125       data Z1,Z2,Z3,Z4 /1600.,1615.,1700.,1715./
 126       data Ncirc /173/              !Half L0 circuits Y>0 (TOT=346 )
 127 c
 128 c
 129 c
 130 c METHOD :  coding from down to up for X
 131 c           coding only Y>0 (if Y<0, NUM-->-NUM)
 132 c
 133 c with NUM = CODE of TRIGGER CIRCUIT
 134 c
 135 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 136 cXXXXXXXXXX
 137 c
 138 c MC1 plane 1
 139 c X11max higher X limit of the circuit
 140 c X11min lower  X limit of the circuit
 141 c
 142 c MC1 plane 2
 143 c X12max=X11max*Z2/Z1
 144 c X12min=X11min*Z2/Z1
 145 c
 146 c MC2 plane 1
 147 c X21max=(X11max+8*XWI_p)*Z3/Z1     ex : CODEX=221 --> XWI_c = 2.12
 148 c X21min=(X11min-8*XWI_m)*Z3/Z1                    --> XWI_p = 1.06
 149 c                                                  --> XWI_m = 2.12
 150 c                                   ex : CODEX=224 --> XWI_c = 2.12
 151 c                                                  --> XWI_p = 4.24
 152 c                                                  --> XWI_m = 2.12
 153 c                                   ex : CODEX=440 --> XWI_c = 4.24
 154 c                                                  --> XWI_p = 0.
 155 c                                                  --> XWI_m = 4.24
 156 c
 157 c MC2 plane 2
 158 c X22max=(X11max+8*XWI_p)*Z4/Z1
 159 c X22min=(X11min-8*XWI_m)*Z4/Z1
 160 c
 161 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 162 cYYYYYYYYYYY
 163 c
 164 c MC1 plane 1
 165 c Y11max  higher Y limit of the circuit
 166 c Y11min  lower  Y limit of the circuit
 167 c XY11max higher X limit of the circuit (in Y)
 168 c XY11min lower  X limit of the circuit (in Y)
 169 c
 170 c MC1 plane 2
 171 c Y12max=Y11max*Z2/Z1
 172 c Y12min=Y11min*Z2/Z1
 173 c XY12max=XY11max*Z2/Z1
 174 c XY12min=XY11min*Z2/Z1
 175 c
 176 c MC2 plane 1
 177 c Y21max=Y11max*Z3/Z1
 178 c Y21min=Y11min*Z3/Z1
 179 c XY21max(n)=XY11max(n)  *Z3/Z1     if IXYCO_p=0
 180 c            XY11max(n+1)*Z3/Z1     if IXYCO_p=1
 181 c XY21min(n)=XY11min(n)  *Z3/Z1     if IXYCO_m=0
 182 c            XY11min(n-1)*Z3/Z1     if IXYCO_m=1
 183 c
 184 c MC2 plane 2
 185 c Y22max=Y11max*Z4/Z1
 186 c Y22min=Y11min*Z4/Z1
 187 c XY22max(n)=XY11max(n)  *Z4/Z1     if IXYCO_p=0
 188 c            XY11max(n+1)*Z4/Z1     if IXYCO_p=1
 189 c XY22min(n)=XY11min(n)  *Z4/Z1     if IXYCO_m=0
 190 c            XY11min(n-1)*Z4/Z1     if IXYCO_m=1
 191 c
 192 c
 193 c CODEY=C1*1000+C2-*100+C2c*10+C2+*1
 194 c      with C1=2/4 for strip width 2.12/4.24 on MC1
 195 c           C2-,C2c,C2+=0/1/2/3 on MC2
 196 c           0--> nothing 1-->4strips 2-->8strips 3-->16strips
 197 c ex : CODEY=2332
 198 c      C1=2.12 (MC1 plane 1)
 199 c      C2-= 16 strips (MC2)
 200 c      C2c= 16 strips (MC2)
 201 c      C2+= 8  strips (MC2)
 202 c      IXYCO_p, IXYCO_m : 0/1 if C2+,C2- =/# 0
 203 c
 204 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 205 c
 206 c   numerotation 1/2 plane (-xxx for other 1/2 plane)
 207 c
 208 c     X
 209 c     ^
 210 c     !
 211 c     ! 121 221 321 421 521 621 715 815 915 1008
 212 c     ! 120                 620 714 814 914 1007
 213 c     ! 119                 619 713 813 913 1007
 214 c     ! 118                 618 712 812 912 1006
 215 c     ! 117                 617 711 811 911 1006
 216 c     ! 116                 616 710 810 910 1005
 217 c     ! 115                 615 710 810 910 1005
 218 c     ! 114 214 314 414     614 709 809 909 1005
 219 c     !         313 413 513 613 709 809 909 1005
 220 c     !             412 512 612 708 808 908 1004
 221 c     ! Beam        411 511 611 708 808 908 1004   --------> Y
 222 c     ! Shield      410 510 610 707 807 907 1004
 223 c     !             409 509 609 707 807 907 1004
 224 c     !         308 408 508 608 706 806 906 1003
 225 c     ! 107 207 307 407     607 706 906 906 1003
 226 c     ! 106                 606 705 805 905 1003
 227 c     ! 105                 605 705 805 905 1003
 228 c     ! 104                 604 704 804 904 1002
 229 c     ! 103                 603 703 803 903 1002
 230 c     ! 102                 602 702 802 902 1001
 231 c     ! 101                 601 701 801 901 1001
 232 c     ! 100 200 300 400 500 600 700 800 900 1000
 233 c
 234 c
 235 cccccccccccccccccccccccccc DATAs CIRCUITS 100-107 114-121 ccccccccccccccccccc
 236
 237        data (NUM(NN),NN=1,16) /100,101,102,103,104,105,106,107,114,115,
 238      +                 116,117,118,119,120,121/
 239
 240        data (CODEX(NN),NN=1,16) /042,422,222,222,221,211,111,110,
 241      +                   011,111,112,122,222,222,224,240/
 242
 243        data (CODEY(NN),NN=1,16) /4011,4110,4012,2122,2222,
 244      +                    2220,2020,2020,
 245      +                   2020,2020,2022,2222,2221,4210,4011,4110/
 246
 247        data (xcmax(1,NN),NN=1,16) /-238,-204,-170,-136,-102,-85,
 248      + -68,-51,
 249      +                     68,85,102,136,170,204,238,306/
 250
 251        data (xcmin(1,NN),NN=1,16) /-306,-238,-204,-170,-136,
 252      +    -102,-85,-68,
 253      +                     51,68,85,102,136,170,204,238/
 254
 255        data (ycmax(1,NN),NN=1,16) /17,17,17,17,17,17,17,17,
 256      +                     17,17,17,17,17,17,17,17/
 257
 258        data (ycmin(1,NN),NN=1,16) /0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0/
 259
 260
 261       data (xycmax(1,NN),NN=1,16) /-238,-170,-170,-136,-102,
 262      +           -51,-51,-51,
 263      +                     102,102,102,136,170,238,238,306/
 264
 265       data (xycmin(1,NN),NN=1,16) /-306,-238,-238,-170,-136,
 266      +                      -102,-102,-102,
 267      +                     51,51,51,102,136,170,170,238/
 268
 269
 270 cccccccccccccccccccccccccc DATAs CIRCUITS 200-207 214-221 ccccccccccccccccccc
 271
 272       data (NUM(NN),NN=17,32) /200,201,202,203,204,205,206,207,214,215,
 273      +                 216,217,218,219,220,221/
 274
 275       data (CODEX(NN),NN=17,32) /042,422,222,222,221,211,111,110,
 276      +                   011,111,112,122,222,222,224,240/
 277
 278       data (CODEY(NN),NN=17,32) /4011,4110,4012,2122,2222,
 279      +      2220,2020,2020,
 280      +                   2020,2020,2022,2222,2221,4210,4011,4110/
 281
 282       data (xcmax(1,NN),NN=17,32) /-238,-204,-170,-136,-102,
 283      +                            -85,-68,-51,
 284      +                    68,85,102,136,170,204,238,306/
 285
 286       data (xcmin(1,NN),NN=17,32) /-306,-238,-204,-170,-136,
 287      +            -102,-85,-68,
 288      +                     51,68,85,102,136,170,204,238/
 289
 290       data (ycmax(1,NN),NN=17,32) /34,34,34,34,34,34,34,
 291      +                    34,34,34,34,34,34,34,34,34/
 292
 293       data (ycmin(1,NN),NN=17,32) /17,17,17,17,17,17,17,
 294      +                    17,17,17,17,17,17,17,17,17/
 295
 296       data (xycmax(1,NN),NN=17,32) /-238,-170,-170,-136,-102,
 297      +              -51,-51,-51,
 298      +                      102,102,102,136,170,238,238,306/
 299
 300       data (xycmin(1,NN),NN=17,32) /-306,-238,-238,-170,-136,
 301      +                     -102,-102,-102,
 302      +                     51,51,51,102,136,170,170,238/
 303
 304
 305 cccccccccccccccccccccccccc DATAs CIRCUITS 300-308 313-321 ccccccccccccccccccc
 306
 307        data (NUM(NN),NN=33,50) /300,301,302,303,304,305,306,307,308,
 308      +                  313,314,315,316,317,318,319,320,321/
 309
 310        data (CODEX(NN),NN=33,50) /042,422,222,222,221,211,111,111,110,
 311      +                    011,111,111,112,122,222,222,224,240/
 312
 313        data (CODEY(NN),NN=33,50) /4011,4110,4012,2122,2222,2220,
 314      + 2022,2220,2020,
 315      +                    2020,2022,2220,2022,2222,2221,4210,4011,4110/
 316
 317        data (xcmax(1,NN),NN=33,50) /-238,-204,-170,-136,-102,-85,
 318      +                   -68,-51,-34,
 319      +                     51,68,85,102,136,170,204,238,306/
 320
 321        data (xcmin(1,NN),NN=33,50) /-306,-238,-204,-170,-136,
 322      +  -102,-85,-68,-51,
 323      +                     34,51,68,85,102,136,170,204,238/
 324
 325        data (ycmax(1,NN),NN=33,50) /51,51,51,51,51,51,51,
 326      +                     51,51,51,51,51,51,51,51,51,51,51/
 327
 328        data (ycmin(1,NN),NN=33,50) /34,34,34,34,34,34,34,
 329      +                     34,34,34,34,34,34,34,34,34,34,34/
 330
 331       data (xycmax(1,NN),NN=33,50) /-238,-170,-170,-136,-102,
 332      +                   -68,-68,-34,-34,
 333      +                     68,68,102,102,136,170,238,238,306/
 334
 335       data (xycmin(1,NN),NN=33,50) /-306,-238,-238,-170,-136,
 336      +                        -102,-102,-68,-68,
 337      +                     34,34,68,68,102,136,170,170,238/
 338
 339
 340 cccccccccccccccccccccccccc DATAs CIRCUITS 400 - 421 ccccccccccccccccccccccccc
 341
 342        data (NUM(NN),NN=51,72) /400,401,402,403,404,405,406,
 343      + 407,408,409,410,
 344      +                  411,412,413,414,415,416,417,418,419,420,421/
 345
 346        data (CODEX(NN),NN=51,72) /042,422,222,222,221,211,111,111,111,
 347      + 111,111,
 348      +                    111,111,111,111,111,112,122,222,222,224,240/
 349
 350        data (CODEY(NN),NN=51,72) /4011,4110,4012,2122,2222,2220,
 351      +  2022,2220,2022,
 352      +                    2220,2022,2220,2022,
 353      +                    2220,2022,2220,2022,2222,2221,4210,4011,4110/
 354
 355        data (xcmax(1,NN),NN=51,72) /-238,-204,-170,-136,-102,-85,
 356      +                 -68,-51,-34,
 357      +                     -17,0,17,34,
 358      +                     51,68,85,102,136,170,204,238,306/
 359
 360        data (xcmin(1,NN),NN=51,72) /-306,-238,-204,-170,-136,-102,
 361      +                    -85,-68,-51,
 362      +                     -34,-17,0,17,
 363      +                     34,51,68,85,102,136,170,204,238/
 364
 365        data (ycmax(1,NN),NN=51,72) /68,68,68,68,68,68,68,68,68,68,68,
 366      +                     68,68,68,68,68,68,68,68,68,68,68/
 367
 368        data (ycmin(1,NN),NN=51,72) /51,51,51,51,51,51,51,51,51,51,51,
 369      +                     51,51,51,51,51,51,51,51,51,51,51/
 370
 371       data (xycmax(1,NN),NN=51,72) /-238,-170,-170,-136,-102,-68,
 372      +        -68,-34,-34,
 373      +                     0,0,34,34,
 374      +                     68,68,102,102,136,170,238,238,306/
 375
 376       data (xycmin(1,NN),NN=51,72) /-306,-238,-238,-170,-136,
 377      +                     -102,-102,-68,-68,
 378      +                     -34,-34,0,0,
 379      +                     34,34,68,68,102,136,170,170,238/
 380
 381
 382 cccccccccccccccccccccccccc DATAs CIRCUITS 500 - 521 ccccccccccccccccccccccccc
 383
 384       data (NUM(NN),NN=73,94) /500,501,502,503,504,505,506,
 385      +                   507,508,509,510,
 386      +                  511,512,513,514,515,516,517,518,519,520,521/
 387
 388       data (CODEX(NN),NN=73,94) /042,422,222,222,221,211,111,
 389      +            111,111,111,111,
 390      +                    111,111,111,111,111,112,122,222,222,224,240/
 391
 392       data (CODEY(NN),NN=73,94) /4011,4110,4012,2122,2222,
 393      +               2220,2022,2220,2022,
 394      +                    2220,2022,2220,2022,
 395      +                    2220,2022,2220,2022,2222,2221,4210,4011,4110/
 396
 397       data (xcmax(1,NN),NN=73,94) /-238,-204,-170,-136,
 398      +                     -102,-85,-68,-51,-34,
 399      +                     -17,0,17,34,
 400      +                     51,68,85,102,136,170,204,238,306/
 401
 402       data (xcmin(1,NN),NN=73,94) /-306,-238,-204,-170,-136,
 403      +                -102,-85,-68,-51,
 404      +                     -34,-17,0,17,
 405      +                     34,51,68,85,102,136,170,204,238/
 406
 407       data (ycmax(1,NN),NN=73,94) /85,85,85,85,85,85,85,85,85,85,85,
 408      +                     85,85,85,85,85,85,85,85,85,85,85/
 409
 410       data (ycmin(1,NN),NN=73,94) /68,68,68,68,68,68,68,68,68,68,68,
 411      +                     68,68,68,68,68,68,68,68,68,68,68/
 412
 413       data (xycmax(1,NN),NN=73,94) /-238,-170,-170,-136,-102,
 414      +            -68,-68,-34,-34,
 415      +                      0,0,34,34,
 416      +                      68,68,102,102,136,170,238,238,306/
 417
 418       data (xycmin(1,NN),NN=73,94) /-306,-238,-238,-170,
 419      +               -136,-102,-102,-68,
 420      +                      -68,-34,-34,0,0,
 421      +                      34,34,68,68,102,136,170,170,238/
 422
 423
 424 cccccccccccccccccccccccccc DATAs CIRCUITS 600 - 621 ccccccccccccccccccccccccc
 425
 426       data (NUM(NN),NN=95,116)/600,601,602,603,604,605,606,
 427      +           607,608,609,610,
 428      +                  611,612,613,614,615,616,617,618,619,620,621/
 429
 430       data (CODEX(NN),NN=95,116)/042,422,222,222,221,211,111,
 431      +      111,111,111,111,
 432      +                    111,111,111,111,111,112,122,222,222,224,240/
 433
 434       data (CODEY(NN),NN=95,116)/4011,4110,4012,2122,2222,
 435      +                    2220,2022,2220,2022,
 436      +                    2220,2022,2220,2022,
 437      +                    2220,2022,2220,2022,2222,2221,4210,4011,4110/
 438
 439       data (xcmax(1,NN),NN=95,116)/-238,-204,-170,-136,
 440      +                     -102,-85,-68,-51,-34,
 441      +                     -17,0,17,34,
 442      +                     51,68,85,102,136,170,204,238,306/
 443
 444       data (xcmin(1,NN),NN=95,116)/-306,-238,-204,-170,-136,
 445      +                    -102,-85,-68,-51,
 446      +                     -34,-17,0,17,
 447      +                     34,51,68,85,102,136,170,204,238/
 448
 449       data (ycmax(1,NN),NN=95,116)/102,102,102,102,102,102,
 450      +               102,102,102,102,
 451      +               102,102,102,102,102,102,102,102,102,102,102,102/
 452
 453       data (ycmin(1,NN),NN=95,116)/85,85,85,85,85,85,85,85,85,85,85,
 454      +                     85,85,85,85,85,85,85,85,85,85,85/
 455
 456       data (xycmax(1,NN),NN=95,116)/-238,-170,-170,-136,-102,-68,
 457      +                -68,-34,-34,
 458      +                      0,0,34,34,
 459      +                      68,68,102,102,136,170,238,238,306/
 460
 461       data (xycmin(1,NN),NN=95,116)/-306,-238,-238,-170,-136,
 462      +              -102,-102,-68,
 463      +                      -68,-34,-34,0,0,
 464      +                      34,34,68,68,102,136,170,170,238/
 465
 466
 467 cccccccccccccccccccccccccc DATAs CIRCUITS 700 - 715 ccccccccccccccccccccccccc
 468
 469       data (NUM(NN),NN=117,132)/700,701,702,703,704,705,706,707,
 470      +                  708,709,710,711,712,713,714,715/
 471
 472       data (CODEX(NN),NN=117,132)/042,422,222,222,222,222,222,222,
 473      +                    222,222,222,222,222,222,224,240/
 474
 475       data (CODEY(NN),NN=117,132)/4022,4220,4023,2233,2333,
 476      +         2333,2333,2333,
 477      +                    2333,2333,2333,2333,2332,4320,4022,4220/
 478
 479       data (xcmax(1,NN),NN=117,132)/-238,-204,-170,-136,-102,-68,-34,0,
 480      +                     34,68,102,136,170,204,238,306/
 481
 482       data (xcmin(1,NN),NN=117,132)/-306,-238,-204,-170,-136,
 483      +  -102,-68,-34,
 484      +                     0,34,68,102,136,170,204,238/
 485
 486       data (ycmax(1,NN),NN=117,132)/136,136,136,136,136,136,136,136,
 487      +                     136,136,136,136,136,136,136,136/
 488
 489       data (ycmin(1,NN),NN=117,132)/102,102,102,102,102,102,102,102,
 490      +                     102,102,102,102,102,102,102,102/
 491
 492       data (xycmax(1,NN),NN=117,132)/-238,-170,-170,
 493      +            -136,-102,-68,-34,0,
 494      +                      34,68,102,136,170,238,238,306/
 495
 496       data (xycmin(1,NN),NN=117,132)/-306,-238,-238,-170,
 497      +            -136,-102,-68,-34,
 498      +                      0,34,68,102,136,170,170,238/
 499
 500
 501
 502 cccccccccccccccccccccccccc DATAs CIRCUITS 800 - 815 ccccccccccccccccccccccccc
 503
 504       data (NUM(NN),NN=133,148)/800,801,802,803,804,805,806,807,
 505      +                  808,809,810,811,812,813,814,815/
 506
 507       data (CODEX(NN),NN=133,148)/042,422,222,222,222,222,222,222,
 508      +                    222,222,222,222,222,222,224,240/
 509
 510       data (CODEY(NN),NN=133,148)/4022,4220,4023,2233,2333,
 511      +              2333,2333,2333,
 512      +                    2333,2333,2333,2333,2332,4320,4022,4220/
 513
 514       data (xcmax(1,NN),NN=133,148)/-238,-204,-170,-136,-102,-68,-34,0,
 515      +                     34,68,102,136,170,204,238,306/
 516
 517       data (xcmin(1,NN),NN=133,148)/-306,-238,-204,-170,
 518      +-136,-102,-68,-34,
 519      +                     0,34,68,102,136,170,204,238/
 520
 521       data (ycmax(1,NN),NN=133,148)/170,170,170,170,170,170,170,170,
 522      +                     170,170,170,170,170,170,170,170/
 523
 524       data (ycmin(1,NN),NN=133,148)/136,136,136,136,136,136,136,136,
 525      +                     136,136,136,136,136,136,136,136/
 526
 527       data (xycmax(1,NN),NN=133,148)/-238,-170,-170,-136,
 528      +-102,-68,-34,0,
 529      +                      34,68,102,136,170,238,238,306/
 530
 531       data (xycmin(1,NN),NN=133,148)/-306,-238,-238,-170,
 532      +      -136,-102,-68,-34,
 533      +                      0,34,68,102,136,170,170,238/
 534
 535
 536 cccccccccccccccccccccccccc DATAs CIRCUITS 900 - 915 ccccccccccccccccccccccccc
 537
 538       data (NUM(NN),NN=149,164)/900,901,902,903,904,905,906,907,
 539      +                  908,909,910,911,912,913,914,915/
 540
 541       data (CODEX(NN),NN=149,164)/042,422,222,222,222,222,222,222,
 542      +                    222,222,222,222,222,222,224,240/
 543
 544       data (CODEY(NN),NN=149,164)/4022,4220,4022,4220,
 545      +4022,4220,4022,4220,
 546      +                    4022,4220,4022,4220,4022,4220,4022,4220/
 547
 548       data (xcmax(1,NN),NN=149,164)/-238,-204,-170,-136,-102,-68,
 549      +-34,0,
 550      +                     34,68,102,136,170,204,238,306/
 551
 552       data (xcmin(1,NN),NN=149,164)/-306,-238,-204,-170,-136,-102,
 553      +              -68,-34,
 554      +                     0,34,68,102,136,170,204,238/
 555
 556       data (ycmax(1,NN),NN=149,164)/204,204,204,204,204,204,204,204,
 557      +                     204,204,204,204,204,204,204,204/
 558
 559       data (ycmin(1,NN),NN=149,164)/170,170,170,170,170,170,170,170,
 560      +                     170,170,170,170,170,170,170,170/
 561
 562       data (xycmax(1,NN),NN=149,164)/-238,-170,-170,-102,-102,
 563      +-34,-34,34,
 564      +                      34,102,102,170,170,238,238,306/
 565
 566       data (xycmin(1,NN),NN=149,164)/-306,-238,-238,-170,-170,
 567      +-102,-102,-34,
 568      +                      -34,34,34,102,102,170,170,238/
 569
 570
 571 cccccccccccccccccccccccccc DATAs CIRCUITS 1000 - 1008 ccccccccccccccccccccccc
 572
 573       data (NUM(NN),NN=165,173)/1000,1001,1002,1003,1004,1005,
 574      +                  1006,1007,1008/
 575
 576       data (CODEX(NN),NN=165,173)/044,444,444,444,444,444,444,444,440/
 577
 578       data (CODEY(NN),NN=165,173)/4033,4333,4333,4333,4333,4333,
 579      +                    4333,4333,4330/
 580
 581       data (xcmax(1,NN),NN=165,173)/-238,-170,-102,-34,
 582      +                     34,102,170,238,306/
 583
 584       data (xcmin(1,NN),NN=165,173)/-306,-238,-170,-102,-34,
 585      +                     34,102,170,238/
 586
 587       data (ycmax(1,NN),NN=165,173)/272,272,272,272,272,272,
 588      +272,272,272/
 589
 590       data (ycmin(1,NN),NN=165,173)/204,204,204,204,204,204,
 591      +204,204,204/
 592
 593       data (xycmax(1,NN),NN=165,173)/-238,-170,-102,-34,
 594      +                      34,102,170,238,306/
 595
 596       data (xycmin(1,NN),NN=165,173)/-306,-238,-170,-102,-34,
 597      +                      34,102,170,238/
 598
 599 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
 600
 601 calculate other quantities associated to a circuit on the other planes
 602 c
 603
 604       do i=1,ncirc                           !loop on all circuits (Y>0)
 605 c
 606 c DECODE CODEX
 607 c
 608       ixwi_m=CODEX(i)/100                    !Dividing integers
 609       ixwi_c=(CODEX(i)-ixwi_m*100)/10
 610       ixwi_p=CODEX(i)-ixwi_m*100-ixwi_c*10
 611 c
 612       if(ixwi_m.eq.0)then
 613        xwi_m(i)=0.
 614       elseif(ixwi_m.eq.1)then
 615        xwi_m(i)=1.0625
 616       elseif(ixwi_m.eq.2)then
 617        xwi_m(i)=2.125
 618       elseif(ixwi_m.eq.4)then
 619        xwi_m(i)=4.25
 620       else
 621        print *,'WARNING :: BAD CODEX value'
 622       endif
 623
 624       if(ixwi_c.eq.0)then
 625        xwi_c(i)=0.
 626       elseif(ixwi_c.eq.1)then
 627        xwi_c(i)=1.0625
 628       elseif(ixwi_c.eq.2)then
 629        xwi_c(i)=2.125
 630       elseif(ixwi_c.eq.4)then
 631        xwi_c(i)=4.25
 632       else
 633        print *,'WARNING :: BAD CODEX value'
 634       endif
 635
 636       if(ixwi_p.eq.0)then
 637        xwi_p(i)=0.
 638       elseif(ixwi_p.eq.1)then
 639        xwi_p(i)=1.0625
 640       elseif(ixwi_p.eq.2)then
 641        xwi_p(i)=2.125
 642       elseif(ixwi_p.eq.4)then
 643        xwi_p(i)=4.25
 644       else
 645        print *,'WARNING :: BAD CODEX value'
 646       endif
 647
 648 c ready to calculate XGEOM parameters
 649
 650       xcmax(2,i)=Z2/Z1*(xcmax(1,i))
 651       xcmin(2,i)=Z2/Z1*(xcmin(1,i))
 652 c
 653 c
 654       xcmax(3,i)=Z3/Z1*(xcmax(1,i)+8*XWI_p(i))
 655       xcmin(3,i)=Z3/Z1*(xcmin(1,i)-8*XWI_m(i))
 656 c
 657       xcmax(4,i)=Z4/Z1*(xcmax(1,i)+8*XWI_p(i))
 658       xcmin(4,i)=Z4/Z1*(xcmin(1,i)-8*XWI_m(i))
 659 c
 660 c
 661 c DECODE CODEY
 662 c
 663       iywi=CODEY(i)/1000                 !Dividing integers
 664       iywi_m=(CODEY(i)-iywi*1000)/100
 665       iywi_c=(CODEY(i)-iywi*1000-iywi_m*100)/10
 666       iywi_p=CODEY(i)-iywi*1000-iywi_m*100-iywi_c*10
 667
 668       if(iywi.eq.2)then
 669        ywi(i)=2.125
 670       elseif(iywi.eq.4)then
 671        ywi(i)=4.25
 672       else
 673        print *,'WARNING :: BAD CODEY value'
 674       endif
 675
 676       if(iywi_m.eq.0)then
 677        ixyco_m=0
 678        nstrip_m(i)=0
 679       elseif(iywi_m.eq.1)then
 680        ixyco_m=1
 681        nstrip_m(i)=4
 682       elseif(iywi_m.eq.2)then
 683        ixyco_m=1
 684        nstrip_m(i)=8
 685       elseif(iywi_m.eq.3)then
 686        ixyco_m=1
 687        nstrip_m(i)=16
 688       else
 689        print *,'WARNING :: BAD CODEY value'
 690       endif
 691
 692       if(iywi_c.eq.0)then
 693        nstrip_c(i)=0
 694       elseif(iywi_c.eq.1)then
 695        nstrip_c(i)=4
 696       elseif(iywi_c.eq.2)then
 697        nstrip_c(i)=8
 698       elseif(iywi_c.eq.3)then
 699        nstrip_c(i)=16
 700       else
 701        print *,'WARNING :: BAD CODEY value'
 702       endif
 703
 704       if(iywi_p.eq.0)then
 705        ixyco_p=0
 706        nstrip_p(i)=0
 707       elseif(iywi_p.eq.1)then
 708        ixyco_p=1
 709        nstrip_p(i)=4
 710       elseif(iywi_p.eq.2)then
 711        ixyco_p=1
 712        nstrip_p(i)=8
 713       elseif(iywi_p.eq.3)then
 714        ixyco_p=1
 715        nstrip_p(i)=16
 716       else
 717        print *,'WARNING :: BAD CODEY value'
 718       endif
 719
 720 c ready to calculate YGEOM parameters
 721
 722
 723       ycmax(2,i)=Z2/Z1*(ycmax(1,i))
 724       ycmin(2,i)=Z2/Z1*(ycmin(1,i))
 725       xycmax(2,i)=Z2/Z1*(xycmax(1,i))
 726       xycmin(2,i)=Z2/Z1*(xycmin(1,i))
 727
 728       ycmax(3,i)=Z3/Z1*(ycmax(1,i))
 729       ycmin(3,i)=Z3/Z1*(ycmin(1,i))
 730
 731       if(ixyco_p.eq.0)then
 732        xycmax(3,i)=Z3/Z1*(xycmax(1,i))
 733       else
 734        xycmax(3,i)=Z3/Z1*(xycmax(1,i+1))
 735       endif
 736
 737       if(ixyco_m.eq.0)then
 738        xycmin(3,i)=Z3/Z1*(xycmin(1,i))
 739       else
 740        xycmin(3,i)=Z3/Z1*(xycmin(1,i-1))
 741       endif
 742
 743       ycmax(4,i)=Z4/Z1*(ycmax(1,i))
 744       ycmin(4,i)=Z4/Z1*(ycmin(1,i))
 745
 746       if(ixyco_p.eq.0)then
 747        xycmax(4,i)=Z4/Z1*(xycmax(1,i))
 748       else
 749        xycmax(4,i)=Z4/Z1*(xycmax(1,i+1))
 750       endif
 751
 752       if(ixyco_m.eq.0)then
 753        xycmin(4,i)=Z4/Z1*(xycmin(1,i))
 754       else
 755        xycmin(4,i)=Z4/Z1*(xycmin(1,i-1))
 756       endif
 757 c
 758 c
 759       enddo               !loop on trigger boards
 760 c
 761       return
 762       end
 763 c
 764 cc
 765       subroutine REMPL(x,y,nhitmax)
 766 c Associate hits to circuits (bit pattern)
 767
 768       real x(4,1000),y(4,1000)  ! 1000=nhitmax
 769
 770       real xu(346),yu(346)                 ! something   in x / y
 771       integer*4 bitpx(4,346),bitpy(4,346)  ! bit-pattern in x / y code DECIMAL
 772       common /hitcirc/xu,yu,bitpx,bitpy
 773
 774 c      integer NUM(173),CODEX(173),CODEY(173)           !173=Ncirc
 775 c      real xcmax(4,173),xcmin(4,173),
 776 c     +     ycmax(4,173),ycmin(4,173),
 777 c     +     xycmax(4,173),xycmin(4,173)
 778 c      real xwi_c(173),xwi_m(173),xwi_p(173),ywi(173)
 779 c      integer nstrip_c(173),nstrip_p(173),nstrip_m(173)
 780 c
 781 c      common/TRIGMAP/xcmax,xcmin,ycmax,ycmin,xycmax,xycmin,
 782 c     +               CODEX,CODEY,
 783 c     +               xwi_c,xwi_m,xwi_p,ywi,
 784 c     +               nstrip_c,nstrip_m,nstrip_p,
 785 c     +               NUM,Z1,Z2,Z3,Z4,Ncirc
 786       common/TRIGMAP/xcmax(4,173),xcmin(4,173),ycmax(4,173),
 787      +               ycmin(4,173),xycmax(4,173),xycmin(4,173),
 788      +               CODEX(173),CODEY(173),
 789      +               xwi_c(173),xwi_m(173),xwi_p(173),ywi(173),
 790      +               nstrip_c(173),nstrip_m(173),nstrip_p(173),
 791      +               NUM(173),Z1,Z2,Z3,Z4,Ncirc
 792       integer nstrip_c,nstrip_p,nstrip_m
 793       integer NUM,CODEX,CODEY,Ncirc
 794
 795
 796       common/debug/idebug
 797
 798       character*32 i2bin
 799
 800
 801 cinit
 802       do icircT=1,(ncirc*2)
 803        xu(icircT)=0
 804        yu(icircT)=0
 805        do iplan=1,4
 806         bitpx(iplan,icircT)=0
 807         bitpy(iplan,icircT)=0
 808        enddo
 809       enddo
 810 c
 811 c fill the hits in histograms
 812 c
 813 c      do iplan=1,4
 814 c       do ihit=1,nhitmax
 815 c        if (x(iplan,ihit).ne.0.or.y(iplan,ihit).ne.0)
 816 c     + call hfill(iplan,x(iplan,ihit),y(iplan,ihit),1.)
 817 c       enddo
 818 c      enddo
 819 c
 820 c fill the circuits in X and Y
 821 c
 822       do iplan=1,4
 823        do ihit=1,nhitmax
 824         if (x(iplan,ihit).ne.0.or.y(iplan,ihit).ne.0)then !save CPU time
 825 c           print *,iplan,ihit,x(iplan,ihit),y(iplan,ihit)
 826 c
 827 ccccccccccccccccccccc BIT PATTERN in X cccccccccccccccccccccccccccc
 828 c
 829         do icirc=1,ncirc
 830          yabs=abs(y(iplan,ihit))
 831          if(x(iplan,ihit).ge.xcmin(iplan,icirc).and.
 832      +      x(iplan,ihit).lt.xcmax(iplan,icirc))then
 833           if(yabs.ge.ycmin(iplan,icirc).and.
 834      +       yabs.lt.ycmax(iplan,icirc))then
 835            if(y(iplan,ihit).ge.0)xu(icirc)=1
 836            if(y(iplan,ihit).lt.0)xu(icirc+ncirc)=1
 837 c
 838            if(iplan.eq.1)then
 839             X11=x(iplan,ihit)-xcmin(iplan,icirc)
 840             nx11=ifix(X11/xwi_c(icirc))+1              !no of hitten strip
 841
 842             if(nX11.le.0.or.nX11.ge.17)
 843      +         print *,'WARNING : NX11 out of range',NX11
 844
 845             if(y(iplan,ihit).ge.0)
 846      +         call sbit1(bitpx(iplan,icirc),nx11)     !sets bit pos nx11 =1
 847             if(y(iplan,ihit).lt.0)
 848      +         call sbit1(bitpx(iplan,icirc+ncirc),nx11)
 849
 850
 851             if (idebug.eq.2)then
 852             print *,' '
 853             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
 854             print *,'x y xcmin xcmax ycmin ycmax ',x(iplan,ihit),
 855      +       y(iplan,ihit),xcmin(iplan,icirc),xcmax(iplan,icirc),
 856      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
 857             print *,'X11 nX11 xwidth=',X11,nX11,xwi_c(icirc)
 858             print *,'bitpx(iplan,icirc)         ',
 859      +        i2bin(bitpx(iplan,icirc),16)
 860             print *,'bitpx(iplan,icirc+ncirc)   ',
 861      +        i2bin(bitpx(iplan,icirc+ncirc),16)
 862             endif
 863
 864            elseif(iplan.eq.2)then
 865 c
 866             X12=x(iplan,ihit)-xcmin(iplan,icirc)
 867             nx12=ifix(X12/(xwi_c(icirc)*Z2/Z1))+1      !no of hitten strip
 868             if(nX12.le.0.or.nX12.ge.17)
 869      +         print *,'WARNING : NX12 out of range',NX12
 870
 871              if(y(iplan,ihit).ge.0)
 872      +          call sbit1(bitpx(iplan,icirc),nx12)    !sets bit pos nx12 =1
 873              if(y(iplan,ihit).lt.0)
 874      +          call sbit1(bitpx(iplan,icirc+ncirc),nx12)
 875
 876
 877             if (idebug.eq.2)then
 878             print *,' '
 879             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
 880             print *,'x y xcmin xcmax ycmin ycmax ',x(iplan,ihit),
 881      +       y(iplan,ihit),xcmin(iplan,icirc),xcmax(iplan,icirc),
 882      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
 883             print *,'X12 nX12 xwidth',X12,nX12,xwi_c(icirc)
 884             print *,'bitpx(iplan,icirc)         ',
 885      +        i2bin(bitpx(iplan,icirc),16)
 886             print *,'bitpx(iplan,icirc+ncirc)   ',
 887      +        i2bin(bitpx(iplan,icirc+ncirc),16)
 888             endif
 889
 890            elseif(iplan.eq.3)then
 891 c
 892             X21=x(iplan,ihit)-xcmin(iplan,icirc)
 893              if(X21.le.(xwi_m(icirc)*8*Z3/Z1))then
 894               nX21=ifix(X21/(xwi_m(icirc)*Z3/Z1))+1    !no of hitten strip
 895              elseif
 896      +        (X21.le.((xwi_m(icirc)*8+xwi_c(icirc)*16)*Z3/Z1))then
 897               XX21=X21-xwi_m(icirc)*8*Z3/Z1
 898               nX21=ifix(XX21/(xwi_c(icirc)*Z3/Z1))+1+8
 899              else
 900               XXX21=X21-(xwi_m(icirc)*8+xwi_c(icirc)*16)*Z3/Z1
 901               nX21=ifix(XXX21/(xwi_p(icirc)*Z3/Z1))+1+8+16
 902              endif
 903
 904             if(nX21.le.0.or.nX21.ge.33)
 905      +         print *,'WARNING : NX21 out of range',NX21
 906
 907             if(y(iplan,ihit).ge.0)
 908      +         call sbit1(bitpx(iplan,icirc),nX21)     !sets bit pos nX21 =1
 909             if(y(iplan,ihit).lt.0)
 910      +         call sbit1(bitpx(iplan,icirc+ncirc),nX21)
 911
 912
 913             if (idebug.eq.2)then
 914             print *,' '
 915             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
 916             print *,'x y xcmin xcmax ycmin ycmax ',x(iplan,ihit),
 917      +       y(iplan,ihit),xcmin(iplan,icirc),xcmax(iplan,icirc),
 918      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
 919             print *,'X21 nX21 xwidth-c xwidth-m xwidth-p =',
 920      +       X21,nX21,xwi_c(icirc),xwi_m(icirc),xwi_p(icirc)
 921             print *,'bitpx(iplan,icirc)         ',
 922      +        i2bin(bitpx(iplan,icirc),32)
 923             print *,'bitpx(iplan,icirc+ncirc)   ',
 924      +        i2bin(bitpx(iplan,icirc+ncirc),32)
 925             endif
 926
 927            elseif(iplan.eq.4)then
 928 c
 929             x22=x(iplan,ihit)-xcmin(iplan,icirc)
 930              if(x22.le.(xwi_m(icirc)*8*Z4/Z1))then
 931               nx22=ifix(x22/(xwi_m(icirc)*Z4/Z1))+1    !no of hitten strip
 932              elseif
 933      +        (x22.le.((xwi_m(icirc)*8+xwi_c(icirc)*16)*Z4/Z1))then
 934               Xx22=x22-xwi_m(icirc)*8*Z4/Z1
 935               nx22=ifix(Xx22/(xwi_c(icirc)*Z4/Z1))+1+8
 936              else
 937               XXx22=x22-(xwi_m(icirc)*8+xwi_c(icirc)*16)*Z4/Z1
 938               nx22=ifix(XXx22/(xwi_p(icirc)*Z4/Z1))+1+8+16
 939              endif
 940
 941             if(nx22.le.0.or.nx22.ge.33)
 942      +         print *,'WARNING : Nx22 out of range',Nx22
 943
 944             if(y(iplan,ihit).ge.0)
 945      +         call sbit1(bitpx(iplan,icirc),nx22)     !sets bit pos nx22 =1
 946             if(y(iplan,ihit).lt.0)
 947      +         call sbit1(bitpx(iplan,icirc+ncirc),nx22)
 948
 949
 950
 951             if (idebug.eq.2)then
 952             print *,' '
 953             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
 954             print *,'x y xcmin xcmax ycmin ycmax ',x(iplan,ihit),
 955      +       y(iplan,ihit),xcmin(iplan,icirc),xcmax(iplan,icirc),
 956      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
 957             print *,'x22 nx22 xwidth-c xwidth-m xwidth-p =',
 958      +       x22,nx22,xwi_c(icirc),xwi_m(icirc),xwi_p(icirc)
 959             print *,'bitpx(iplan,icirc)         ',
 960      +        i2bin(bitpx(iplan,icirc),32)
 961             print *,'bitpx(iplan,icirc+ncirc)   ',
 962      +        i2bin(bitpx(iplan,icirc+ncirc),32)
 963             endif
 964
 965
 966            endif
 967           endif
 968          endif
 969
 970 ccccccccccccccccccccc BIT PATTERN in Y cccccccccccccccccccccccccccc
 971
 972          if(x(iplan,ihit).ge.xycmin(iplan,icirc).and.
 973      +      x(iplan,ihit).lt.xycmax(iplan,icirc))then
 974           if(yabs.ge.ycmin(iplan,icirc).and.
 975      +       yabs.lt.ycmax(iplan,icirc))then
 976            if(y(iplan,ihit).ge.0)yu(icirc)=1
 977            if(y(iplan,ihit).lt.0)yu(icirc+ncirc)=1
 978 c
 979            if(iplan.eq.1)then
 980             Y11=yabs-ycmin(iplan,icirc)
 981             nY11=ifix(Y11/ywi(icirc))+1                !no of hitten strip
 982             nstrip=nstrip_c(icirc)
 983             if(nY11.le.0.or.nY11.ge.(nstrip+1))
 984      +         print *,'WARNING : NY11 out of range',NY11
 985
 986             if(y(iplan,ihit).ge.0)
 987      +         call sbit1(bitpy(iplan,icirc),nY11)     !sets bit pos nY11 =1
 988             if(y(iplan,ihit).lt.0)
 989      +         call sbit1(bitpy(iplan,icirc+ncirc),nY11)
 990
 991
 992             if (idebug.eq.2)then
 993             print *,' '
 994             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
 995             print *,'x y xycmin xycmax ycmin ycmax ',x(iplan,ihit),
 996      +       y(iplan,ihit),xycmin(iplan,icirc),xycmax(iplan,icirc),
 997      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
 998             print *,'Y11 nY11 nstripy ywidth=',
 999      +               Y11,nY11,nstrip,ywi(icirc)
1000             print *,'bitpy(iplan,icirc)         ',
1001      +        i2bin(bitpy(iplan,icirc),nstrip)
1002             print *,'bitpy(iplan,icirc+ncirc)   ',
1003      +        i2bin(bitpy(iplan,icirc+ncirc),nstrip)
1004              endif
1005
1006
1007            elseif(iplan.eq.2)then
1008 c
1009             Y12=yabs-ycmin(iplan,icirc)
1010             nY12=ifix(Y12/(ywi(icirc)*Z2/Z1))+1        !no of hitten strip
1011             nstrip=nstrip_c(icirc)
1012             if(nY12.le.0.or.nY12.ge.(nstrip+1))
1013      +         print *,'WARNING : NY12 out of range',NY12
1014
1015             if(y(iplan,ihit).ge.0)
1016      +         call sbit1(bitpy(iplan,icirc),nY12)     !sets bit pos nY12 =1
1017             if(y(iplan,ihit).lt.0)
1018      +         call sbit1(bitpy(iplan,icirc+ncirc),nY12)
1019
1020
1021             if (idebug.eq.2)then
1022             print *,' '
1023             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
1024             print *,'x y xycmin xycmax ycmin ycmax ',x(iplan,ihit),
1025      +       y(iplan,ihit),xycmin(iplan,icirc),xycmax(iplan,icirc),
1026      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
1027             print *,'Y12 nY12 nstripy ywidth=',
1028      +               Y12,nY12,nstrip,ywi(icirc)
1029             print *,'bitpy(iplan,icirc)         ',
1030      +        i2bin(bitpy(iplan,icirc),nstrip)
1031             print *,'bitpy(iplan,icirc+ncirc)   ',
1032      +        i2bin(bitpy(iplan,icirc+ncirc),nstrip)
1033             endif
1034
1035            elseif(iplan.eq.3)then
1036 c
1037             Y21=yabs-ycmin(iplan,icirc)
1038             nstrip=nstrip_c(icirc)
1039             nY21=ifix(Y21/(ywi(icirc)*Z3/Z1))+1        !no of hitten strip
1040              if(nY21.le.0.or.nY21.ge.(nstrip+1))
1041      +         print *,'WARNING : NY21 out of range',NY21
1042
1043              if(y(iplan,ihit).ge.0)
1044      +         call sbit1(bitpy(iplan,icirc),nY21)     !sets bit pos nY21 =1
1045              if(y(iplan,ihit).lt.0)
1046      +         call sbit1(bitpy(iplan,icirc+ncirc),nY21)
1047
1048 c various Y strip width in the same circuit
1049             if(nstrip_p(icirc).ne.0.and.nstrip_p(icirc).lt.nstrip)then
1050              xyextrap=xycmax(1,icirc)*Z3/Z1
1051               if(x(iplan,ihit).gt.xyextrap)then
1052                if((mod(ny21,2)).eq.1)then
1053                 if(y(iplan,ihit).ge.0)
1054      +           call sbit1(bitpy(iplan,icirc),nY21+1) !sets bit pos nY21+1 =1
1055                 if(y(iplan,ihit).lt.0)
1056      +           call sbit1(bitpy(iplan,icirc+ncirc),nY21+1)
1057                endif
1058                if((mod(ny21,2)).eq.0)then
1059                 if(y(iplan,ihit).ge.0)
1060      +           call sbit1(bitpy(iplan,icirc),nY21-1) !sets bit pos nY21-1 =1
1061                 if(y(iplan,ihit).lt.0)
1062      +           call sbit1(bitpy(iplan,icirc+ncirc),nY21-1)
1063                endif
1064               endif
1065             endif
1066
1067             if(nstrip_m(icirc).ne.0.and.nstrip_m(icirc).lt.nstrip)then
1068              xyextram=xycmin(1,icirc)*Z3/Z1
1069               if(x(iplan,ihit).lt.xyextram)then
1070                if((mod(ny21,2)).eq.1)then
1071                 if(y(iplan,ihit).ge.0)
1072      +           call sbit1(bitpy(iplan,icirc),nY21+1) !sets bit pos nY21+1 =1
1073                 if(y(iplan,ihit).lt.0)
1074      +           call sbit1(bitpy(iplan,icirc+ncirc),nY21+1)
1075                endif
1076                if((mod(ny21,2)).eq.0)then
1077                 if(y(iplan,ihit).ge.0)
1078      +           call sbit1(bitpy(iplan,icirc),nY21-1) !sets bit pos nY21-1 =1
1079                 if(y(iplan,ihit).lt.0)
1080      +           call sbit1(bitpy(iplan,icirc+ncirc),nY21-1)
1081                endif
1082               endif
1083             endif
1084
1085             if (idebug.eq.2)then
1086             print *,' '
1087             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
1088             print *,'x y xycmin xycmax ycmin ycmax ',x(iplan,ihit),
1089      +       y(iplan,ihit),xycmin(iplan,icirc),xycmax(iplan,icirc),
1090      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
1091             print *,'Y21 nY21 nstripc ywidth= nstripm nstripp',
1092      +      Y21,nY21,nstrip,ywi(icirc),nstrip_m(icirc),nstrip_p(icirc)
1093             print *,'bitpy(iplan,icirc)         ',
1094      +        i2bin(bitpy(iplan,icirc),nstrip)
1095             print *,'bitpy(iplan,icirc+ncirc)   ',
1096      +        i2bin(bitpy(iplan,icirc+ncirc),nstrip)
1097             if(nstrip_p(icirc).ne.0.and.nstrip_p(icirc).lt.nstrip)then
1098             print *,'special changt de largeur Y : p plus larges'
1099             print *,'xyextrap, mod(ny21,2)=',xyextrap, mod(ny21,2)
1100             endif
1101             if(nstrip_m(icirc).ne.0.and.nstrip_m(icirc).lt.nstrip)then
1102             print *,'special changt de largeur Y : m plus larges'
1103             print *,'xyextram, mod(ny21,2)=',xyextrap, mod(ny21,2)
1104             endif
1105             endif
1106
1107            elseif(iplan.eq.4)then
1108 c
1109             Y22=yabs-ycmin(iplan,icirc)
1110             nstrip=nstrip_c(icirc)
1111              nY22=ifix(Y22/(ywi(icirc)*Z4/Z1))+1       !no of hitten strip
1112              if(nY22.le.0.or.nY22.ge.(nstrip+1))
1113      +         print *,'WARNING : NY22 out of range',NY22
1114
1115              if(y(iplan,ihit).ge.0)
1116      +         call sbit1(bitpy(iplan,icirc),nY22)     !sets bit pos nY22 =1
1117              if(y(iplan,ihit).lt.0)
1118      +         call sbit1(bitpy(iplan,icirc+ncirc),nY22)
1119 c
1120 c various Y strip width in the same circuit
1121             if(nstrip_p(icirc).ne.0.and.nstrip_p(icirc).lt.nstrip)then
1122              xyextrap=xycmax(1,icirc)*Z4/Z1
1123               if(x(iplan,ihit).gt.xyextrap)then
1124                if((mod(ny22,2)).eq.1)then
1125                 if(y(iplan,ihit).ge.0)
1126      +           call sbit1(bitpy(iplan,icirc),nY22+1) !sets bit pos nY22+1 =1
1127                 if(y(iplan,ihit).lt.0)
1128      +           call sbit1(bitpy(iplan,icirc+ncirc),nY22+1)
1129                endif
1130                if((mod(ny22,2)).eq.0)then
1131                 if(y(iplan,ihit).ge.0)
1132      +           call sbit1(bitpy(iplan,icirc),nY22-1) !sets bit pos nY22-1 =1
1133                 if(y(iplan,ihit).lt.0)
1134      +           call sbit1(bitpy(iplan,icirc+ncirc),nY22-1)
1135                endif
1136               endif
1137             endif
1138
1139             if(nstrip_m(icirc).ne.0.and.nstrip_m(icirc).lt.nstrip)then
1140              xyextram=xycmin(1,icirc)*Z4/Z1
1141               if(x(iplan,ihit).lt.xyextram)then
1142                if((mod(ny22,2)).eq.1)then
1143                 if(y(iplan,ihit).ge.0)
1144      +           call sbit1(bitpy(iplan,icirc),nY22+1) !sets bit pos nY22+1 =1
1145                 if(y(iplan,ihit).lt.0)
1146      +           call sbit1(bitpy(iplan,icirc+ncirc),nY22+1)
1147                endif
1148                if((mod(ny22,2)).eq.0)then
1149                 if(y(iplan,ihit).ge.0)
1150      +           call sbit1(bitpy(iplan,icirc),nY22-1) !sets bit pos nY22-1 =1
1151                 if(y(iplan,ihit).lt.0)
1152      +           call sbit1(bitpy(iplan,icirc+ncirc),nY22-1)
1153                endif
1154               endif
1155             endif
1156
1157
1158             if (idebug.eq.2)then
1159             print *,' '
1160             print *,'iplan,ihit,icirc,NUM=',iplan,ihit,icirc,num(icirc)
1161             print *,'x y xycmin xycmax ycmin ycmax ',x(iplan,ihit),
1162      +       y(iplan,ihit),xycmin(iplan,icirc),xycmax(iplan,icirc),
1163      +       ycmin(iplan,icirc),ycmax(iplan,icirc)
1164             print *,'Y22 nY22 nstripc ywidth= nstripm nstripp',
1165      +      Y22,nY22,nstrip,ywi(icirc),nstrip_m(icirc),nstrip_p(icirc)
1166             print *,'bitpy(iplan,icirc)         ',
1167      +        i2bin(bitpy(iplan,icirc),nstrip)
1168             print *,'bitpy(iplan,icirc+ncirc)   ',
1169      +        i2bin(bitpy(iplan,icirc+ncirc),nstrip)
1170             if(nstrip_p(icirc).ne.0.and.nstrip_p(icirc).lt.nstrip)then
1171             print *,'special changt de largeur Y : p plus larges'
1172             print *,'xyextrap, mod(ny22,2)=',xyextrap,mod(ny22,2)
1173             endif
1174             if(nstrip_m(icirc).ne.0.and.nstrip_m(icirc).lt.nstrip)then
1175             print *,'special changt de largeur Y : m plus larges'
1176             print *,'xyextram, mod(ny22,2)=',xyextrap, mod(ny22,2)
1177             endif
1178             endif
1179
1180
1181            endif                                       !iplan=1-4
1182 c
1183           endif
1184          endif
1185 c
1186         enddo                                          !circuit
1187         endif                                          !x-y ne 0
1188        enddo                                           !hits
1189       enddo                                            !plan
1190
1191
1192       return
1193       end
1194 c
1195 cc
1196       subroutine ALGO(itrigR,itrigL0,itrigH0)
1197
1198
1199 cc
1200 c      integer NUM(173),CODEX(173),CODEY(173)           !173=Ncirc
1201 c      real xcmax(4,173),xcmin(4,173),
1202 c     +     ycmax(4,173),ycmin(4,173),
1203 c     +     xycmax(4,173),xycmin(4,173)
1204 c      real xwi_c(173),xwi_m(173),xwi_p(173),ywi(173)
1205 c      integer nstrip_c(173),nstrip_p(173),nstrip_m(173)
1206 c
1207 c      common/TRIGMAP/xcmax,xcmin,ycmax,ycmin,xycmax,xycmin,
1208 c     +               CODEX,CODEY,
1209 c     +               xwi_c,xwi_m,xwi_p,ywi,
1210 c     +               nstrip_c,nstrip_m,nstrip_p,
1211 c     +               NUM,Z1,Z2,Z3,Z4,Ncirc
1212       common/TRIGMAP/xcmax(4,173),xcmin(4,173),ycmax(4,173),
1213      +               ycmin(4,173),xycmax(4,173),xycmin(4,173),
1214      +               CODEX(173),CODEY(173),
1215      +               xwi_c(173),xwi_m(173),xwi_p(173),ywi(173),
1216      +               nstrip_c(173),nstrip_m(173),nstrip_p(173),
1217      +               NUM(173),Z1,Z2,Z3,Z4,Ncirc
1218       integer nstrip_c,nstrip_p,nstrip_m
1219       integer NUM,CODEX,CODEY,Ncirc
1220
1221       common/debug/idebug
1222
1223 cc
1224       real xu(346),yu(346)                 ! something   in x / y
1225       integer*4 bitpx(4,346),bitpy(4,346)  ! bit-pattern in x / y code DECIMAL
1226
1227       common /hitcirc/xu,yu,bitpx,bitpy
1228
1229       character*32 i2bin
1230
1231       integer*4 dble1_0,dble1_1_14,dble1_15
1232       integer*4 dble2_0,dble2_1_30,dble2_31
1233       integer*4 sgle1_0,sgle1_1_14,sgle1_15
1234       integer*4 sgle2_0,sgle2_1_30,sgle2_31
1235       integer*4 sgle1A,sgle1B,sgle2A,sgle2B
1236       integer*4 x1_1,x1_2a,x1_2b,x1_2c,x2_1,x2_2a,x2_2b,x2_2c
1237       integer*4 dblex1(346),dblex2(346)
1238       integer*4 sglex1(346),sglex2(346)
1239       integer*4 dbley1(346),dbley2(346)
1240       integer*4 sgley1(346),sgley2(346)
1241       integer*4 thrl(346)
1242       integer*4 co_l(16,346),co_y(16,346)
1243       integer*4 ib,jb,iabit,ibbit,icbit,ibit_l,ibit_y
1244       integer*4 sign_l(346),val_y(346)
1245       integer*4 sign_lv
1246 c L0 L.U.T. (et L2)
1247       integer*4 dnp,dnm,dsup,dinf
1248       integer*4 dev_2(346),num_x2(346),num_y2(346),signdev_2(346)
1249       integer*4 devmin(16),stripnum(16),idevmin(16)
1250       real yL2(346),x1L2(346),x2L2(346)
1251 c
1252       JBIT (IZW,IZP)     = IBITS (IZW,IZP-1,1)
1253       JBYT (IZW,IZP,NZB) = IBITS (IZW,IZP-1,NZB)
1254
1255
1256
1257 c Threshold Road MAX +/-8
1258       data thrl /346*131071/  ! +/- 8 strips
1259 c      data thrl /346*4064/  ! +/- 3 strips
1260 c
1261 c Datas for L.U.T. calculations
1262       data zF /975./
1263       data ptcalLow,ptcalHigh /.60,1.6/
1264 c
1265 c init
1266 c
1267       do k=1,346
1268          sign_l(k)=0
1269          val_y(k)=0
1270          yL2(k)=0
1271          x1L2(k)=0
1272          x2L2(k)=0
1273       enddo
1274 c
1275 c histogram number and distribution of hitten circuits
1276 c
1277
1278       nHxcirc=0
1279       nHycirc=0
1280       nHcirc=0
1281       do icircT=1,ncirc*2
1282        if(xu(icircT).eq.1)nHxcirc=nHxcirc+1
1283        if(yu(icircT).eq.1)nHycirc=nHycirc+1
1284        if(xu(icircT).eq.1.or.yu(icircT).eq.1)nHcirc=nHcirc+1  !circuit non-vide
1285        if(xu(icircT).eq.1.or.yu(icircT).eq.1)then
1286         if(icircT.le.ncirc)then
1287 c         call hf1(11,float(icircT),1.)
1288         else
1289 c         call hf1(11,-float(icircT-ncirc),1.)
1290         endif
1291        endif
1292       enddo
1293 c      call hf1(10,float(nHcirc),1.)
1294 c      call hf1(12,float(nHxcirc),1.)
1295 c      call hf1(13,float(nHycirc),1.)
1296
1297 c
1298 c Calculate Singles and Doubles (including MINI-ROADS)
1299 c
1300       do icircT=1,ncirc*2
1301       if(xu(icircT).eq.1)then  !loop on non-empty X circuit
1302 c
1303 c memo : integer=ibits(integer,start_pos,lenght), start_pos_min=0
1304 c
1305 c MC1
1306        x1_1=ibits(bitpx(1,icircT),0,1)
1307        x1_2a=ibits(bitpx(2,icircT),0,1)
1308        x1_2b=ibits(bitpx(2,icircT),1,1)
1309        dble1_0=iand(x1_1,ior(x1_2a,x1_2b))
1310
1311        x1_1 =ibits(bitpx(1,icircT),1,14)
1312        x1_2a=ibits(bitpx(2,icircT),0,14)
1313        x1_2b=ibits(bitpx(2,icircT),1,14)
1314        x1_2c=ibits(bitpx(2,icircT),2,14)
1315        dble1_1_14=iand(x1_1,ior(x1_2a,ior(x1_2b,x1_2c)))
1316
1317        x1_1=ibits(bitpx(1,icircT),15,1)
1318        x1_2a=ibits(bitpx(2,icircT),14,1)
1319        x1_2b=ibits(bitpx(2,icircT),15,1)
1320        dble1_15=iand(x1_1,ior(x1_2a,x1_2b))
1321
1322        dblex1(icircT)=dble1_0+2*dble1_1_14+(2**15)*dble1_15  !dble X on MC1
1323
1324        x1_1=ibits(bitpx(1,icircT),0,16)
1325        sgle1A=ieor(x1_1,dblex1(icircT))
1326
1327        x1_2a=ibits(dblex1(icircT),0,1)
1328        x1_2b=ibits(dblex1(icircT),1,1)
1329        sgle1_0=ior(x1_2a,x1_2b)
1330
1331        x1_2a=ibits(dblex1(icircT),0,14)
1332        x1_2b=ibits(dblex1(icircT),1,14)
1333        x1_2c=ibits(dblex1(icircT),2,14)
1334        sgle1_1_14=ior(x1_2a,ior(x1_2b,x1_2c))
1335
1336        x1_2a=ibits(dblex1(icircT),14,1)
1337        x1_2b=ibits(dblex1(icircT),15,1)
1338        sgle1_15=ior(x1_2a,x1_2b)
1339
1340        sgle1B=sgle1_0+2*sgle1_1_14+(2**15)*sgle1_15
1341        sgle1B=not(sgle1B)
1342        sgle1B=iand(sgle1B,bitpx(2,icircT))
1343
1344        sglex1(icircT)=ior(sgle1A,sgle1B)                     !sgle X on MC1
1345
1346
1347
1348        dbley1(icircT)=iand(bitpy(1,icircT),bitpy(2,icircT))     !dble Y on MC1
1349        sgley1(icircT)=
1350      + ieor(ior(bitpy(1,icircT),bitpy(2,icircT)),dbley1(icircT))!sgle Y on MC1
1351
1352
1353        if(idebug.eq.3)then
1354          if(icircT.le.ncirc)then
1355           nocirc=NUM(icircT)
1356          else
1357           nocirc=-NUM(icircT-ncirc)
1358          endif
1359          print *,' '
1360          print *,'no circuit ',icircT,nocirc
1361          print *,'bitpx(1,icircT)         ',i2bin(bitpx(1,icircT),16)
1362          print *,'bitpx(2,icircT)         ',i2bin(bitpx(2,icircT),16)
1363          print *,'dblex1(icircT)          ',i2bin(dblex1(icircT),16)
1364          print *,'sglex1(icircT)          ',i2bin(sglex1(icircT),16)
1365          print *,' '
1366          print *,'bitpy(1,icircT)         ',i2bin(bitpy(1,icircT),16)
1367          print *,'bitpy(2,icircT)         ',i2bin(bitpy(2,icircT),16)
1368          print *,'dbley1(icircT)          ',i2bin(dbley1(icircT),16)
1369          print *,'sgley1(icircT)          ',i2bin(sgley1(icircT),16)
1370
1371        endif
1372
1373
1374 c MC2
1375        x2_1=ibits(bitpx(3,icircT),0,1)
1376        x2_2a=ibits(bitpx(4,icircT),0,1)
1377        x2_2b=ibits(bitpx(4,icircT),1,1)
1378        dble2_0=iand(x2_1,ior(x2_2a,x2_2b))
1379
1380        x2_1 =ibits(bitpx(3,icircT),1,30)
1381        x2_2a=ibits(bitpx(4,icircT),0,30)
1382        x2_2b=ibits(bitpx(4,icircT),1,30)
1383        x2_2c=ibits(bitpx(4,icircT),2,30)
1384        dble2_1_30=iand(x2_1,ior(x2_2a,ior(x2_2b,x2_2c)))
1385
1386        x2_1=ibits(bitpx(3,icircT),31,1)
1387        x2_2a=ibits(bitpx(4,icircT),30,1)
1388        x2_2b=ibits(bitpx(4,icircT),31,1)
1389        dble2_31=iand(x2_1,ior(x2_2a,x2_2b))
1390
1391        dblex2(icircT)=dble2_0+2*dble2_1_30+(2.**31)*dble2_31  !dble X on MC2
1392
1393        x2_1=ibits(bitpx(3,icircT),0,32)
1394        sgle2A=ieor(x2_1,dblex2(icircT))
1395
1396        x2_2a=ibits(dblex2(icircT),0,1)
1397        x2_2b=ibits(dblex2(icircT),1,1)
1398        sgle2_0=ior(x2_2a,x2_2b)
1399
1400        x2_2a=ibits(dblex2(icircT),0,30)
1401        x2_2b=ibits(dblex2(icircT),1,30)
1402        x2_2c=ibits(dblex2(icircT),2,30)
1403        sgle2_1_30=ior(x2_2a,ior(x2_2b,x2_2c))
1404
1405        x2_2a=ibits(dblex2(icircT),30,1)
1406        x2_2b=ibits(dblex2(icircT),31,1)
1407        sgle2_31=ior(x2_2a,x2_2b)
1408
1409        sgle2B=sgle2_0+2*sgle2_1_30+(2.**31)*sgle2_31
1410        sgle2B=not(sgle2B)
1411        sgle2B=iand(sgle2B,bitpx(4,icircT))
1412
1413        sglex2(icircT)=ior(sgle2A,sgle2B)                      !sgle X on MC2
1414
1415
1416        dbley2(icircT)=iand(bitpy(3,icircT),bitpy(4,icircT))     !dble Y on MC2
1417        sgley2(icircT)=
1418      + ieor(ior(bitpy(3,icircT),bitpy(4,icircT)),dbley2(icircT))!sgle Y on MC2
1419
1420
1421        if(idebug.eq.3)then
1422          print*,' '
1423          print *,'bitpx(3,icircT) ',i2bin(bitpx(3,icircT),32)
1424          print *,'bitpx(4,icircT) ',i2bin(bitpx(4,icircT),32)
1425          print *,'dblex2(icircT)  ',i2bin(dblex2(icircT),32)
1426          print *,'sglex2(icircT)  ',i2bin(sglex2(icircT),32)
1427          print *,' '
1428          print *,'bitpy(3,icircT) ',i2bin(bitpy(3,icircT),16)
1429          print *,'bitpy(4,icircT) ',i2bin(bitpy(4,icircT),16)
1430          print *,'dbley2(icircT)  ',i2bin(dbley2(icircT),16)
1431          print *,'sgley2(icircT)  ',i2bin(sgley2(icircT),16)
1432        endif
1433
1434        endif                    !non-empty X circuit
1435       enddo                     !circuits
1436 c
1437 c
1438 c method DS (X only) : reduce the sensibility of the trigger to single hits
1439 c                      without loosing signal
1440 c
1441        do icircT=1,ncirc*2
1442         if(xu(icircT).eq.1)then  !loop on non-empty X circuit to save CPU time
1443          if(dblex1(icircT).ne.0)sglex1(icircT)=0
1444          if(dblex2(icircT).ne.0)sglex2(icircT)=0
1445         endif
1446        enddo
1447 c
1448 c
1449 c
1450 c coincidence 3/4 in +/- 8 strips road
1451 c
1452       do icircT=1,2*ncirc
1453       if(xu(icircT).eq.1)then  !loop on non-empty X circuit
1454        do ib=1,16
1455         do jb=ib,ib+16
1456          iabit=iand(jbit(sglex1(icircT),ib),jbit(dblex2(icircT),jb))
1457          ibbit=ior(jbit(sglex2(icircT),jb),jbit(dblex2(icircT),jb))
1458          icbit=ior(iabit,iand(ibbit,jbit(dblex1(icircT),ib)))
1459          ibit_l=iand(icbit,jbit(thrl(icircT),jb-ib+1))
1460          call sbit(ibit_l,co_l(ib,icircT),jb-ib+1)
1461         enddo
1462        enddo
1463
1464        if(idebug.eq.4)then
1465         if (xu(icircT).ne.0)then
1466          if(icircT.le.ncirc)then
1467           nocirc=NUM(icircT)
1468          else
1469           nocirc=-NUM(icircT-ncirc)
1470          endif
1471          print *,' '
1472          print *,'no circuit ',icircT,nocirc
1473          print *,'sglex1(icircT)          ',i2bin(sglex1(icircT),16)
1474          print *,'dblex1(icircT)          ',i2bin(dblex1(icircT),16)
1475          print *,'sglex2(icircT)  ',i2bin(sglex2(icircT),32)
1476          print *,'dblex2(icircT)  ',i2bin(dblex2(icircT),32)
1477          print *,' '
1478          print *,'thrl(icircT)             ',i2bin(thrl(icircT),17)
1479           do ib=1,16
1480          print *,'co_l(ib,icircT)          ',i2bin(co_l(ib,icircT),17)
1481           enddo
1482
1483         endif
1484        endif
1485 c
1486 c signe info (X) before L.U.T.
1487 c
1488        ibit_l=0
1489        do ib=1,16
1490         ibit_l=ior(ibit_l,co_l(ib,icircT))
1491        enddo
1492         ibgauche=ibits(ibit_l,9,8)
1493         ibcentre=ibits(ibit_l,8,1)
1494         ibdroite=ibits(ibit_l,0,8)
1495        if(ibgauche.eq.0.and.ibcentre.eq.0.and.ibdroite.eq.0)then
1496         sign_l(icircT)=0                                !00
1497        elseif(ibgauche.eq.0.and.ibcentre.eq.0)then
1498         sign_l(icircT)=1                                !01
1499        elseif(ibdroite.eq.0.and.ibcentre.eq.0)then
1500         sign_l(icircT)=2                                !10
1501        else
1502         sign_l(icircT)=3                                !11
1503        endif
1504
1505
1506        if(idebug.eq.4)then
1507         if (xu(icircT).ne.0)then
1508          print*,' '
1509          print *,'sign_l(icircT)          ',i2bin(sign_l(icircT),2)
1510         endif
1511        endif
1512       endif                       !non-empty X circuits
1513       enddo                       !circuits
1514 c
1515 c coincidences 3/4 Y (constant ROAD +/- 1)
1516 c
1517       do icircT=1,ncirc*2
1518       if(xu(icircT).eq.1)then  !loop on non-empty X circuit
1519
1520 c
1521          if(icircT.le.ncirc)then
1522           nstripy=nstrip_c(icircT)
1523          else
1524           nstripy=nstrip_c(icircT-ncirc)
1525          endif
1526 c
1527        do ib=1,nstripy
1528         do jb=ib,ib+2                      !+/- 1 strip
1529
1530          iabit=iand(jbit(sgley1(icircT),ib),jbit(2*dbley2(icircT),jb))
1531          ibbit=ior(jbit(2*sgley2(icircT),jb),jbit(2*dbley2(icircT),jb))
1532          ibit_y=ior(iabit,iand(ibbit,jbit(dbley1(icircT),ib)))
1533          call sbit(ibit_y,co_y(ib,icircT),jb-ib+1)
1534         enddo
1535        enddo
1536 c validation Y (OUTPUT L0 Y)
1537        ibit_y=0
1538        do ib=1,nstripy
1539         ibit_y=ior(ibit_y,co_y(ib,icircT))
1540        enddo
1541
1542
1543        if(ibit_y.eq.0)then
1544         val_y(icircT)=0           !nothing in Y
1545        else
1546         val_y(icircT)=1           !something in Y
1547        endif
1548
1549        if(idebug.eq.5)then
1550         if (xu(icircT).ne.0)then
1551          if(icircT.le.ncirc)then
1552           nocirc=NUM(icircT)
1553          else
1554           nocirc=-NUM(icircT-ncirc)
1555          endif
1556          print *,' '
1557          print *,'no circuit ',icircT,nocirc
1558          print *,'sgley1(icircT)  ',i2bin(sgley1(icircT),nstripy)
1559          print *,'dbley1(icircT)  ',i2bin(dbley1(icircT),nstripy)
1560          print *,'sgley2(icircT)  ',i2bin(sgley2(icircT),nstripy)
1561          print *,'dbley2(icircT)  ',i2bin(dbley2(icircT),nstripy)
1562          print *,' '
1563           do ib=1,nstripy
1564          print *,'co_y(ib,icircT)          ',i2bin(co_y(ib,icircT),3)
1565           enddo
1566          print *,' '
1567          print *,'val_y(icircT)          ',i2bin(val_y(icircT),1)
1568
1569         endif
1570        endif
1571       endif                     !non-empty X circuits
1572       enddo                     !circuits
1573 c
1574 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
1575 c infos in +/- 8 road are available : sign_l(346) and val_y(346)
1576 c where sign_l = 2 bits integer (00 nothing, 01 right dev.,
1577 c                                10 left dev., 11 zero dev. or ambiguity)
1578 c and val_y=0 (nothing) or 1 (validation of X)
1579 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
1580 c
1581 c Calculate trigger response in +/-8 strip road (for debugging)
1582 c
1583       igl_0=0
1584       idl_0=0
1585       igh_0=0
1586       idh_0=0
1587       il_0=0
1588       ih_0=0
1589       itrigR=0
1590
1591
1592
1593       do icircT=1,2*ncirc
1594       if(xu(icircT).eq.1)then  !loop on non-empty X circuit
1595
1596 c with Y validation
1597        sign_lv=sign_l(icircT)*val_y(icircT)
1598        if (sign_lv.eq.1.or.sign_lv.eq.3)idl_0=idl_0+1
1599        if (sign_lv.eq.2.or.sign_lv.eq.3)igl_0=igl_0+1
1600        if (sign_lv.ge.1)il_0=il_0+1
1601
1602       endif                   !non-empty X circuits
1603       enddo                   !circuits
1604
1605
1606       if(il_0.ge.2.and.
1607      + idl_0.ne.0.and.igl_0.ne.0)itrigR=1       !U.S. trigger Road +/-8
1608 *       if(il_0.ge.1)itrigR=1
1609 c
1610 c
1611 c
1612 c
1613 c INTERFACE L0-L.U.T. (and L2) in X :
1614 c       nO CIRCUIT + Hitten STRIP (MC1) + DEVMIN (MC2) + SIGNDEV
1615 c
1616 c  L.U.T. (and L2) calculation ARE done only if a LOCAL trigger
1617 c                in +/-8 strips (sign_lv.ne.0) is found
1618 c
1619 c
1620       do icircT=1,ncirc*2
1621       sign_lv=sign_l(icircT)*val_y(icircT)
1622        if(sign_lv.ne.0)then                        !local trigger +/- 8 strips
1623 c
1624         do ib=1,16
1625          if(co_l(ib,icircT).eq.0)then
1626           devmin(ib)=100                            !dummy
1627          elseif(iand(256,co_l(ib,icircT)).ne.0)then
1628           devmin(ib)=0                              !central strip
1629          elseif(iand(128,co_l(ib,icircT)).ne.0)then
1630           devmin(ib)=-1                             ! DOWN priority and <0
1631          elseif(iand(512,co_l(ib,icircT)).ne.0)then
1632           devmin(ib)=1                              ! UP >0 (convention x>0)
1633          elseif(iand(64,co_l(ib,icircT)).ne.0)then
1634           devmin(ib)=-2
1635          elseif(iand(1024,co_l(ib,icircT)).ne.0)then
1636           devmin(ib)=2
1637          elseif(iand(32,co_l(ib,icircT)).ne.0)then
1638           devmin(ib)=-3
1639          elseif(iand(2048,co_l(ib,icircT)).ne.0)then
1640           devmin(ib)=3
1641          elseif(iand(16,co_l(ib,icircT)).ne.0)then
1642           devmin(ib)=-4
1643          elseif(iand(4096,co_l(ib,icircT)).ne.0)then
1644           devmin(ib)=4
1645          elseif(iand(8,co_l(ib,icircT)).ne.0)then
1646           devmin(ib)=-5
1647          elseif(iand(8192,co_l(ib,icircT)).ne.0)then
1648           devmin(ib)=5
1649          elseif(iand(4,co_l(ib,icircT)).ne.0)then
1650           devmin(ib)=-6
1651          elseif(iand(16384,co_l(ib,icircT)).ne.0)then
1652           devmin(ib)=6
1653          elseif(iand(2,co_l(ib,icircT)).ne.0)then
1654           devmin(ib)=-7
1655          elseif(iand(32768,co_l(ib,icircT)).ne.0)then
1656           devmin(ib)=7
1657          elseif(iand(1,co_l(ib,icircT)).ne.0)then
1658           devmin(ib)=-8
1659          elseif(iand(65536,co_l(ib,icircT)).ne.0)then
1660           devmin(ib)=8
1661          endif
1662         enddo
1663 c order deviation (cf CERNLIB)
1664         do ib =1,16
1665           idevmin(ib)=iabs(devmin(ib))
1666         enddo
1667         call SORTZV(idevmin,stripnum,16,-1,0,0)
1668 c
1669         num_x2(icircT)=stripnum(1)
1670         dev_2(icircT)=devmin(stripnum(1))
1671          if(dev_2(icirct).lt.0)then
1672            signdev_2(icircT)=1
1673          elseif(dev_2(icirct).eq.0)then
1674            signdev_2(icircT)=3
1675          else
1676            signdev_2(icircT)=2
1677          endif
1678
1679
1680 c
1681 c
1682 c INTERFACE L0-L.U.T. (and L2) in Y : nO CIRCUIT + Hitten STRIP (MC1)
1683 c
1684          if(icircT.le.ncirc)then
1685           nstripy=nstrip_c(icircT)
1686          else
1687           nstripy=nstrip_c(icircT-ncirc)
1688          endif
1689 c
1690        do ib=nstripy,1,-1       !keep lower weight (convention)
1691         if(co_y(ib,icircT).ne.0)num_y2(icircT)=ib
1692        enddo
1693
1694
1695
1696         if(idebug.eq.6)then
1697           if(icircT.le.ncirc)then
1698            nocirc=NUM(icircT)
1699           else
1700            nocirc=-NUM(icircT-ncirc)
1701           endif
1702           print *,' '
1703           print *,'no circuit ',icircT,nocirc
1704          do ib=1,16
1705          print *,'strip no ',ib,':      ',i2bin(co_l(ib,icircT),17),
1706      +   '      dev_min=',devmin(ib)
1707          enddo
1708          print *,'num_x2(icircT)= ',num_x2(icircT)
1709          print *,' dev_2(icircT)= ',dev_2(icircT)
1710 c
1711          print *,' '
1712          do ib=1,nstripy
1713          print *,'strip no ',ib,':      ',
1714      +   i2bin(co_y(ib,icircT),3)
1715          enddo
1716          print *,'num_y2(icircT)= ',num_y2(icircT)
1717         endif
1718
1719
1720        endif                               !end local trigger +/- 8 strips
1721       enddo
1722
1723 c
1724 c
1725 cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
1726 c Datas for L0===>L.U.T. (and L2) :
1727 c                       icircT
1728 c                       num_x2(icircT)===>X MC1
1729 c                       dev_2(icircT) ===>X MC2
1730 c                    signdev_2(icircT)===>sign dev_2(icircT)
1731 c                       num_y2(icircT)===>Y MC1
1732 ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
1733 c
1734 c Now calculate X MC1 and X MC2 and Y MC1
1735 c ===> Back to trigger GEOMETRY
1736 c
1737       do icircT=1,2*ncirc
1738       sign_lv=sign_l(icircT)*val_y(icircT)
1739        if(sign_lv.ne.0)then                !local trigger +/- 8 strips
1740 calculation of YL2 (middle of the strip)
1741         if(icircT.le.ncirc)then            !Y>0
1742          yL2(icircT)=ycmin(1,icircT)+ywi(icircT)*(num_y2(icircT)-.5)
1743         else                               ! Y<0
1744          yL2(icircT)=-ycmin(1,icircT-ncirc)
1745      +               -ywi(icircT-ncirc)*(num_y2(icircT)-.5)
1746         endif
1747
1748 calculation of X1L2 (MC1) (middle of the strip)
1749         if(icircT.le.ncirc)then            ! Y>0
1750          x1L2(icircT)=xcmin(1,icircT)+xwi_c(icircT)*(num_x2(icircT)-.5)
1751         else                               ! Y<0
1752          x1L2(icircT)=xcmin(1,icircT-ncirc)
1753      +               +xwi_c(icircT-ncirc)*(num_x2(icircT)-.5)
1754         endif
1755
1756 calculation of X2L2 (MC2) (middle of the strip)
1757
1758           if(icircT.le.ncirc)then
1759            xxwi_c=xwi_c(icircT)
1760            xxwi_p=xwi_p(icircT)
1761            xxwi_m=xwi_m(icircT)
1762           else
1763            xxwi_c=xwi_c(icircT-ncirc)
1764            xxwi_p=xwi_p(icircT-ncirc)
1765            xxwi_m=xwi_m(icircT-ncirc)
1766           endif
1767
1768          if(dev_2(icircT).ge.0)then
1769           Dnp=16-num_x2(icircT)
1770           Dnm=0                                         !RAZ, not used
1771           dsup=dev_2(icircT)-Dnp
1772           dinf=0                                        !RAZ, not used
1773            if(dsup.le.0)then
1774             x2L2(icircT)=(x1L2(icircT)+xxwi_c*dev_2(icircT))*Z3/Z1
1775            else
1776             x2L2(icircT)=(x1L2(icircT)+xxwi_c*Dnp+xxwi_p*dsup)*Z3/Z1
1777            endif
1778          else
1779           Dnp=0                                         !RAZ, not used
1780           Dnm=num_x2(icircT)-1
1781           dinf=iabs(dev_2(icircT))-Dnm
1782           dsup=0                                        !RAZ, not used
1783            if(dinf.le.0)then
1784             x2L2(icircT)=(x1L2(icircT)-xxwi_c*iabs(dev_2(icircT)))*Z3/Z1
1785            else
1786             x2L2(icircT)=(x1L2(icircT)-xxwi_c*Dnm-xxwi_m*dinf)*Z3/Z1
1787            endif
1788          endif
1789
1790          if(idebug.eq.7)then
1791           if(icircT.le.ncirc)then
1792            nocirc=NUM(icircT)
1793            yywi=ywi(icircT)
1794           else
1795            nocirc=-NUM(icircT-ncirc)
1796            yywi=ywi(icircT-ncirc)
1797           endif
1798           print *,' '
1799           print *,'no circuit ',icircT,nocirc
1800           print *,'num_y2 ywi yL2  = ',num_y2(icircT),yywi,yL2(icircT)
1801           print *,'num_x2 xwiC x1L2 = ',
1802      +             num_x2(icircT),xxwi_c,x1L2(icircT)
1803           print *,'dev_2 DnP Dnm dsup dinf ',
1804      +             dev_2(icircT),Dnp,Dnm,dsup,dinf
1805           print *,'xwiM xwiC xwiP x2L2 = ',
1806      +             xxwi_m,xxwi_c,xxwi_p,x2L2(icircT)
1807           print *,'sign dev =',i2bin(sign_lv,2)
1808          endif
1809
1810        endif                                  !end local trigger +/- 8 strips
1811       enddo                                   !circuits
1812
1813 c
1814 c
1815 ccc   NOW pt calculation (for the L.U.T)
1816 c
1817 cinit
1818       idlc_2p=0
1819       idhc_2p=0
1820       idlc_2m=0
1821       idhc_2m=0
1822       idlc_2=0
1823       idhc_2=0
1824 c
1825       do icircT=1,2*ncirc
1826        sign_lv=sign_l(icircT)*val_y(icircT)
1827        if(sign_lv.ne.0)then                    !local pt trigger +/- 8 strips
1828 c calculate tetadev (x-z plane)
1829           anum=z3*(x2L2(icircT)-x1L2(icircT))/(z3-z1)
1830           anum=anum-x2L2(icircT)
1831           tetadev=anum/zF
1832 c calculate ptcalc (need rF)
1833           xF=(x2L2(icircT)-x1L2(icircT))*(z3-zF)/(z3-z1)
1834           xF=x2L2(icircT)-xF
1835           yF=yL2(icircT)*zF/z1
1836           rF=sqrt(xF*xF+yF*yF)
1837           if(abs(tetadev).ge.0.00001)then
1838            ptcalc=(rF/zF)*abs(0.9/tetadev)
1839           else                                             ! infinity
1840            ptcalc=10000.
1841           endif
1842
1843 c  for L0 pt cut using the sign of the deviation signdev_2(icircT)
1844
1845        if (ptcalc.ge.ptcalLow.and.
1846      +    (signdev_2(icircT).eq.1.or.signdev_2(icircT).eq.3))
1847      +     idlc_2p=idlc_2p+1
1848        if (ptcalc.ge.ptcalLow.and.
1849      +    (signdev_2(icircT).eq.2.or.signdev_2(icircT).eq.3))
1850      +     idlc_2m=idlc_2m+1
1851        if (ptcalc.ge.ptcalLow.and.signdev_2(icircT).ge.1)
1852      +     idlc_2=idlc_2+1
1853
1854        if (ptcalc.ge.ptcalHigh.and.
1855      +    (signdev_2(icircT).eq.1.or.signdev_2(icircT).eq.3))
1856      +     idhc_2p=idhc_2p+1
1857        if (ptcalc.ge.ptcalHigh.and.
1858      +    (signdev_2(icircT).eq.2.or.signdev_2(icircT).eq.3))
1859      +     idhc_2m=idhc_2m+1
1860        if (ptcalc.ge.ptcalHigh.and.signdev_2(icircT).ge.1)
1861      +     idhc_2=idhc_2+1
1862
1863          if(idebug.eq.8)then
1864           if(icircT.le.ncirc)then
1865            nocirc=NUM(icircT)
1866           else
1867            nocirc=-NUM(icircT-ncirc)
1868           endif
1869           print *,' '
1870           print *,'no circuit ',icircT,nocirc
1871           print *,'compteurs trigger L0 (boucle circuit ) '
1872           print *,'signdev_2 =',i2bin(signdev_2(icircT),2)
1873           print *,'ptcalc,ptcalLow',ptcalc,ptcalLow
1874           print *,'LOW  pt calc + et -, mult:',idlc_2p,idlc_2m,idlc_2
1875           print *,'ptcalc,ptcalHigh',ptcalc,ptcalHigh
1876           print *,'HIGH pt calc + et -, mult:',idhc_2p,idhc_2m,idhc_2
1877          endif
1878
1879
1880        endif                             !end local trigger +/- 8 strips
1881       enddo                              !circuits
1882
1883 c
1884 c for L0 L.U.T. trigger response
1885 c
1886        itrigL0=0
1887        itrigH0=0
1888
1889       if(idlc_2.ge.2.and.
1890      + idlc_2p.ne.0.and.idlc_2m.ne.0)itrigL0=1
1891       if(idhc_2.ge.2.and.
1892      + idhc_2p.ne.0.and.idhc_2m.ne.0)itrigH0=1
1893 *       if(idlc_2.ge.1)itrigL0=1
1894 *       if(idhc_2.ge.1)itrigH0=1
1895 c
1896       return
1897       end
1898 c
1899 c
1900 cc
1901       subroutine seqtest
1902
1903       real x(4,1000),y(4,1000)     ! 1000=nhitmax
1904       common /hitRPC/ x,y
1905 c
1906 c generation of test sequence of hits for debuging proposes
1907 c
1908        x(1,1)= 58.8
1909        y(1,1)=-52.3
1910        x(2,1)= 60.6
1911        y(2,1)=-52.7
1912        x(3,1)= 64.8
1913        y(3,1)=-55.5
1914        x(4,1)= 65.5
1915        y(4,1)=-56.0
1916
1917        x(1,2)=-105.4
1918        y(1,2)=-26.7
1919        x(2,2)=-106.7
1920        y(2,2)=-27.0
1921        x(3,2)=-114.4
1922        y(3,2)=-28.9
1923        x(4,2)=-115.7
1924        y(4,2)=-29.2
1925
1926       return
1927       end
1928 c
1929 c
1930 c
1931 c function that returns a string corresp. to the integer coded in bits
1932 c
1933       CHARACTER*(*) FUNCTION i2bin(iw,kbit)
1934 c      i2bin='' ???????????????
1935
1936       i2bin=' '
1937       do ib= 0,kbit-1
1938        if(ibits(iw,ib,1).eq.0)then
1939         i2bin='0'//i2bin
1940        else
1941         i2bin='1'//i2bin
1942        endif
1943       enddo
1944       return
1945       end
1946
1947       SUBROUTINE SORTZV (A,INDEX,N1,MODE,NWAY,NSORT)
1948 C
1949 C CERN PROGLIB# M101    SORTZV          .VERSION KERNFOR  3.15  820113
1950 C ORIG. 02/10/75
1951 C
1952       DIMENSION A(N1),INDEX(N1)
1953 C
1954 C
1955       N = N1
1956       IF (N.LE.0)            RETURN
1957       IF (NSORT.NE.0) GO TO 2
1958       DO 1 I=1,N
1959     1 INDEX(I)=I
1960 C
1961     2 IF (N.EQ.1)            RETURN
1962       IF (MODE)    10,20,30
1963    10 CALL SORTTI (A,INDEX,N)
1964       GO TO 40
1965 C
1966    20 CALL SORTTC(A,INDEX,N)
1967       GO TO 40
1968 C
1969    30 CALL SORTTF (A,INDEX,N)
1970 C
1971    40 IF (NWAY.EQ.0) GO TO 50
1972       N2 = N/2
1973       DO 41 I=1,N2
1974       ISWAP = INDEX(I)
1975       K = N+1-I
1976       INDEX(I) = INDEX(K)
1977    41 INDEX(K) = ISWAP
1978    50 RETURN
1979       END
1980 *     ========================================
1981       SUBROUTINE SORTTF (A,INDEX,N1)
1982 C
1983       DIMENSION A(N1),INDEX(N1)
1984 C
1985       N = N1
1986       DO 3 I1=2,N
1987       I3 = I1
1988       I33 = INDEX(I3)
1989       AI = A(I33)
1990     1 I2 = I3/2
1991       IF (I2) 3,3,2
1992     2 I22 = INDEX(I2)
1993       IF (AI.LE.A (I22)) GO TO 3
1994       INDEX (I3) = I22
1995       I3 = I2
1996       GO TO 1
1997     3 INDEX (I3) = I33
1998     4 I3 = INDEX (N)
1999       INDEX (N) = INDEX (1)
2000       AI = A(I3)
2001       N = N-1
2002       IF (N-1) 12,12,5
2003     5 I1 = 1
2004     6 I2 = I1 + I1
2005       IF (I2.LE.N) I22= INDEX(I2)
2006       IF (I2-N) 7,9,11
2007     7 I222 = INDEX (I2+1)
2008       IF (A(I22)-A(I222)) 8,9,9
2009     8 I2 = I2+1
2010       I22 = I222
2011     9 IF (AI-A(I22)) 10,11,11
2012    10 INDEX(I1) = I22
2013       I1 = I2
2014       GO TO 6
2015    11 INDEX (I1) = I3
2016       GO TO 4
2017    12 INDEX (1) = I3
2018       RETURN
2019       END
2020 *     ========================================
2021       SUBROUTINE SORTTI (A,INDEX,N1)
2022 C
2023       INTEGER A,AI
2024       DIMENSION A(N1),INDEX(N1)
2025 C
2026       N = N1
2027       DO 3 I1=2,N
2028       I3 = I1
2029       I33 = INDEX(I3)
2030       AI = A(I33)
2031     1 I2 = I3/2
2032       IF (I2) 3,3,2
2033     2 I22 = INDEX(I2)
2034       IF (AI.LE.A (I22)) GO TO 3
2035       INDEX (I3) = I22
2036       I3 = I2
2037       GO TO 1
2038     3 INDEX (I3) = I33
2039     4 I3 = INDEX (N)
2040       INDEX (N) = INDEX (1)
2041       AI = A(I3)
2042       N = N-1
2043       IF (N-1) 12,12,5
2044     5 I1 = 1
2045     6 I2 = I1 + I1
2046       IF (I2.LE.N) I22= INDEX(I2)
2047       IF (I2-N) 7,9,11
2048     7 I222 = INDEX (I2+1)
2049       IF (A(I22)-A(I222)) 8,9,9
2050     8 I2 = I2+1
2051       I22 = I222
2052     9 IF (AI-A(I22)) 10,11,11
2053    10 INDEX(I1) = I22
2054       I1 = I2
2055       GO TO 6
2056    11 INDEX (I1) = I3
2057       GO TO 4
2058    12 INDEX (1) = I3
2059       RETURN
2060       END
2061 *     ========================================
2062       SUBROUTINE SORTTC (A,INDEX,N1)
2063 C
2064       INTEGER A,AI
2065       DIMENSION A(N1),INDEX(N1)
2066 C
2067       N = N1
2068       DO 3 I1=2,N
2069       I3 = I1
2070       I33 = INDEX(I3)
2071       AI = A(I33)
2072     1 I2 = I3/2
2073       IF (I2) 3,3,2
2074     2 I22 = INDEX(I2)
2075       IF(ICMPCH(AI,A(I22)))3,3,21
2076    21 INDEX (I3) = I22
2077       I3 = I2
2078       GO TO 1
2079     3 INDEX (I3) = I33
2080     4 I3 = INDEX (N)
2081       INDEX (N) = INDEX (1)
2082       AI = A(I3)
2083       N = N-1
2084       IF (N-1) 12,12,5
2085     5 I1 = 1
2086     6 I2 = I1 + I1
2087       IF (I2.LE.N) I22= INDEX(I2)
2088       IF (I2-N) 7,9,11
2089     7 I222 = INDEX (I2+1)
2090       IF (ICMPCH(A(I22),A(I222))) 8,9,9
2091     8 I2 = I2+1
2092       I22 = I222
2093     9 IF (ICMPCH(AI,A(I22))) 10,11,11
2094    10 INDEX(I1) = I22
2095       I1 = I2
2096       GO TO 6
2097    11 INDEX (I1) = I3
2098       GO TO 4
2099    12 INDEX (1) = I3
2100       RETURN
2101       END
2102 *     ========================================
2103       FUNCTION ICMPCH(IC1,IC2)
2104 C     FUNCTION TO COMPARE TWO 4 CHARACTER EBCDIC STRINGS - IC1,IC2
2105 C     ICMPCH=-1 IF HEX VALUE OF IC1 IS LESS THAN IC2
2106 C     ICMPCH=0  IF HEX VALUES OF IC1 AND IC2 ARE THE SAME
2107 C     ICMPCH=+1 IF HEX VALUES OF IC1 IS GREATER THAN IC2
2108       I1=IC1
2109       I2=IC2
2110       IF(I1.GE.0.AND.I2.GE.0)GOTO 40
2111       IF(I1.GE.0)GOTO 60
2112       IF(I2.GE.0)GOTO 80
2113       I1=-I1
2114       I2=-I2
2115       IF(I1-I2)80,70,60
2116  40   IF(I1-I2)60,70,80
2117  60   ICMPCH=-1
2118       RETURN
2119  70   ICMPCH=0
2120       RETURN
2121  80   ICMPCH=1
2122       RETURN
2123       END
2124
2125
2126       SUBROUTINE SBIT (IT,IZW,IZP)
2127 C
2128 C CERN PROGLIB# M421    SBIT            .VERSION KERNFOR  4.23  891215
2129 C MOD. true default 24/2/89, JZ
2130 C
2131 C     This non-ANSI code is a default which may be slow, if so
2132 C     it should be replaced by a machine-specific fast routine
2133
2134       IF(IAND(IT,1).EQ.0) THEN
2135          IZW = IBCLR(IZW,IZP-1)
2136       ELSE
2137          IZW = IBSET(IZW,IZP-1)
2138       ENDIF
2139
2140       END
2141
2142       SUBROUTINE SBIT1 (IZW,IZP)
2143 C
2144 C CERN PROGLIB# M421    SBIT1           .VERSION KERNFOR  4.23  891215
2145 C MOD. true default 24/2/89, JZ
2146 C
2147 C     This non-ANSI code is a default which may be slow, if so
2148 C     it should be replaced by a machine-specific fast routine
2149
2150       IZW = IBSET (IZW, IZP-1)
2151
2152       END
2153