PROGRAM ATXT252
C  CALCUL DISTR ANGLES BEETW  PARTICLES
C  FOR EVENT S-32 FOR N MEASUREMENTS
C  fl atxt252.for g77 -o ee atxt252.f  
C  ./ee < ini77.txt > xxx252.res
C  DATA OF LAST REDACTION: 19 fevr 2010 LEPEKHIN
C
      COMMON /EVN/IN(8),JN(4),KI(5,6)
      COMMON /ELEF/F(2000),MM
      COMMON /ELE/XXN(2000),NX
      COMMON /NHI/NHI(1000),IHI(10),NLM,GR(10),SA(10)
      COMMON /AXXX/NNHI,AXX(10),
     *XX(10),XX2(10),X,KM
      COMMON /ISKNUL/ ISKNUL(100),ISK,NNN
      IDAY=19
      MANS=2
      JEAR=2010
      MM=0
      ISK=0
      NX=0
      NNHI=10
      I=1
      NLM=84
      DO 49 K=1,NNHI
      IHI(K)=I
      I=I+NLM
      AXX(K)=0.0
      XX(K)=0.0
      XX2(K)=0.0
      GR(K)=-1000.0
      SA(K)=25.0
   49 CONTINUE 
      SA(1)=25.0
      GR(1)=-1000.0     
      GR(2)=0.0
      SA(2)=25.0
      GR(4)=0.0
      SA(4)=25.0
      GR(5)=0.0
      SA(5)=25.0
      GR(6)=-1000.0
      SA(6)=25.0
      GR(7)=0.0
      SA(7)=25.0
      GR(10)=0.0
      SA(10)=25.0
      GR(8)=-1000.0
      SA(8)=25.0
      GR(9)=0.0
      SA(9)=25.0
      DO 549 K=1,1000
  549 NHI(K)=0
      PRINT 102
      PRINT 112,IDAY,MANS,JEAR
  102 FORMAT (2X,'THIS EVENTS CONTRNTS ONLY ALPHA PARTICLES'/
     *2X,'FROM S-32 200  GEV PER NUCLEON '/
     *2X,'ALL DATA IN ini77.txt '/)
  112 FORMAT (5('*'),2X,I2,'.',I2,'.',I4,2X,'atxt252.f',5('*'))
      PRINT 113
  113 FORMAT (2X,'LIST OF ANGLES > 1000 MICRORAD ISKLUCHENY'/
     *2X,'N IN LIST N V JOURN ANGLES')
C 62
   27 CALL EVN32
      IF (IN(1).EQ.0) GO TO 66
      CALL USR32
      GO TO 27
   66 IF (ISK.EQ.0) GO TO 28
      PRINT 123
  123 FORMAT (2X,'NOMERA SOB PO JURNALU GDE FINOL > 1000 MIKRORAD')
      PRINT 124,(ISKNUL(K),K=1,ISK)
  124 FORMAT (10I5) 
   28 CALL LITPRI
C      DO 13 K=1,NX
C   13 F(K)=XXN(K)
C      MM=NX
C      CALL RANGIR
C      CALL RAVNOM
C      PRINT 121,NX
C  121 FORMAT (2X,'ALL CHISEL=',2X,I8/)
C      PRINT 333,(F(K),K=1,NX)
C  333 FORMAT (F10.3)
      STOP
      END
C   84
      SUBROUTINE USR32
      COMMON /ELE/XXN(2000),NX
      COMMON /EVN/ IN(8),JN(4),KI(5,6)
      COMMON /NHI/ IHI(1000),IHS(10),NLM,GR(10),SA(10)
      COMMON /AXXX/ NNHI,AXX(10),XX(10),XX2(10),X,KM
      COMMON /ISKNUL/ ISKL(100),ISK,NNN
      DIMENSION SRY(5),EE(10)
      DATA CX,CY,Q/1.25,0.02122,1.4142/
      NNN=0
      NFR=IN(4)
      NT=IN(8)
      E=NT
      AY=JN(1)*CY
      AX=JN(2)*CX
      KM=1
      X=(AY/AX)*1000000.0
      IF (ABS(X).LT.1000.0) GO TO 7
      IF (ISK.GT.100) GO TO 7
      NNN=2
      ISK=ISK+1
      ISKL(ISK)=IN(2)
    7 CONTINUE
      CALL HISTO
      PY=X
      KK=NFR+1
      SRX=KI(1,KK)*CX
      DO 4 K=1,NFR
      S=0.0  
      DO 2 L=1,NT
    2 S=S+KI(L,K)*CY
    4 SRY(K)=S/E
      M=NFR-1
      KM=2
      DO 12 K=1,M
      L=K+1
      A=SRY(K)
      DO 11 J=L,NFR
      B=SRY(J)
      X=(ABS(A-B)*1000000.0)/SRX
      X=X/Q
      IF (X.GT.1000.0) GO TO 11
      IF (X.LT.0.0001) GO TO 11
      CALL HISTO
C      GO TO 11
C   13 PRINT 101,IN(1),IN(2),X
C  101 FORMAT (2X,I4,1X,I4,1X,F10.2,1X,'TET I-K') 
   11 CONTINUE
   12 CONTINUE
      IF (NNN.GT.1) GO TO 75
      PY=0.0
      DO 1 K=1,NFR
      X=(SRY(K)/SRX)*1000000.0-PY
      IF (ABS(X).LT.1000.0) GO TO 6
      PRINT 103,IN(1),IN(2),X
  103 FORMAT (2X,I4,1X,I4,1X,F10.2,1X,'TET WITH PERV')
      GO TO 1
    6 CONTINUE
      NX=NX+1
      XXN(NX)=ABS(X)
      KM=3
      CALL HISTO
      X=ABS(X)
      KM=4
      CALL HISTO
    1 CONTINUE
   75 DO 22 L=1,NT
      DO 21 K=1,NFR
      S=KI(L,K)*CY
      KM=6
      X=(S/SRX)*1000000.0
      SRY(K)=S
      IF (ABS(X).GT.1000.0) GO TO 21
      CALL HISTO
      KM=7
      X=ABS(X)
      CALL HISTO
   21 CONTINUE
      M=NFR-1
      KM=5
      DO 24 KK=1,M
      LL=KK+1
      A=SRY(KK)
      DO 25 J=LL,NFR
      B=SRY(J)
      X=(ABS(A-B)*1000000.0)/(SRX*Q)
      IF (X.GT.1000.0) GO TO 25
      CALL HISTO
   25 CONTINUE
   24 CONTINUE
   22 CONTINUE
      S1=250.0
      S2=10.0
      DO 31 K=1,NFR
      CALL GAUSS(S1,A)
      CALL GAUSS(S2,B)
      X=A+B
      KM=8
      CALL HISTO
      EE(K)=ABS(X)/Q
      KM=9
      X=ABS(X)
   31 CALL HISTO
      M=NFR-1
      KM=10
      DO 32 KK=1,M
      LL=KK+1
      A=EE(KK)
      DO 33 J=LL,NFR
      X=ABS(A-EE(J))
   33 CALL HISTO
   32 CONTINUE
      RETURN
      END
C
      SUBROUTINE RAVNOM
      COMMON /ELEF/X(2000),NX
      DIMENSION NRAV(10),DRAV(10),III(10)
      DATA DRAV/31.2,62.5,96.0,130.0,168.0,
     *210.0,260.0,320.0,410.0,750.0/
C      DATA DRAV/17.0,37.5,57.3,73.6,93.6,117.0,128.0,
C     *160.0,195.0,210.0,238.0,270.0,304.0,348.0,
C     *410.0,512.0,748.0/ 
      DO 1 K=1,10
    1 NRAV(K)=0
      DO 31 K=1,NX
      A=X(K)
      DO 30 L=1,10
      G=DRAV(L)
      IF (A.LT.G) NRAV(L)=NRAV(L)+1
   30 CONTINUE
   31 CONTINUE
      NN=NRAV(1)
      III(1)=NN
      DO 41 K=2,10
      III(K)=NRAV(K)-NRAV(K-1)
   41 NN=III(K)+NN
      PRINT 101,NN
  101 FORMAT (2X,I6)
      PRINT 100,(III(K),K=1,10)
  100 FORMAT (2X,'RAVNOMERNOE RASPREDELENIE'/10I5/10I5)
      RETURN
      END
C 108
      SUBROUTINE EVN32
      COMMON /EVN/IN(8),JN(4),KI(5,6)
      READ (5,100) (IN(K),K=1,8)
      N=0
      DO 2 K=1,8
    2 N=N+IN(K)
      IF (N.EQ.0) GO TO 77
      IF (IN(1).GT.1000) GO TO 78
      IF (IN(2).GT.1000) GO TO 78
      IF (IN(3).GT.25) GO TO 78
      IF (IN(4).GT.5) GO TO 78
      IF (IN(5).GT.31) GO TO 78
      IF (IN(6).GT.12) GO TO 78
      IF (IN(7).GT.2010) GO TO 78
      IF (IN(8).GT.5) GO TO 78 
      NT=IN(8)
      NF=IN(4)
      KK=NF+1
      READ (5,100) (JN(K),K=1,4)
      IF (ABS(JN(1)).GT.300) GO TO 78
      IF (JN(2).GT.30000) GO TO 78
      IF (JN(3).GT.200) GO TO 78
      IF (JN(4).GT.200) GO TO 78
  100 FORMAT (8I6)
    4 DO 1 K=1,NT
    1 READ (5,100) (KI(K,L),L=1,KK)
      GO TO 77
   78 PRINT 100,(IN(K),K=1,8)
      PRINT 100,(JN(K),K=1,4)
      GO TO 4   
   77 RETURN
      END
C 122
      SUBROUTINE RANGIR
      COMMON /ELEF/F(2000),MM
      N=MM
      M=N-1
    7 DO 31 L=1,M
      LI=L+1
      A=F(LI)
      IF (A.GT.F(L)) GO TO 31
      F(LI)=F(L)
      F(L)=A
   31 CONTINUE
      M=M-1
      IF (M.GT.0) GO TO 7
      RETURN
      END
C 138
      SUBROUTINE BXXX
      COMMON /ELEF/F(2000),MM
      COMMON /AXXX/NNHI,AXX(10),
     *XX(10),XX2(10),X,I
      XX(I)=XX(I)+X
      XX2(I)=XX2(I)+X*X
      AXX(I)=AXX(I)+1.0
      IF (I.NE.3) GO TO 77
      MM=MM+1
      F(MM)=X
   77 RETURN
      END
C 151
      SUBROUTINE HISTO
      COMMON /NHI/NHI(1000),IHI(10),NLM,GR(10),SA(10)
      COMMON /AXXX/NNHI,AXX(10),
     *XX(10),XX2(10),T,KM
      CALL BXXX
      M=(T-GR(KM))/SA(KM)
      IFL=0
      X=GR(KM)+SA(KM)*(NLM-4)
      IF (T.LT.GR(KM)) GO TO 1
      IF (T.GT.X) GO TO 2
      GO TO 3
    1 M=NLM-4
      IFL=2
      GO TO 3
    2 M=NLM-3
      IFL=2
    3 K=M+IHI(KM)
      NHI(K)=NHI(K)+1
      K=NLM+IHI(KM)-1
      NHI(K)=NHI(K)+1
      K=K-1
      IF (IFL.GT.1) GO TO 77
      NHI(K)=NHI(K)+1
   77 RETURN
      END
C  177
      SUBROUTINE LITPRI
      COMMON /ELEF/F(2000),MM
      COMMON /NHI/IHI(1000),IHS(10),NLM,GR(10),SA(10)
      COMMON /AXXX/NNHI,AXX(10),
     *XX(10),XX2(10),X,KM
      DIMENSION SR(10),SS(10),EE(80)
C 184
      DO 1 K=1,NNHI
C 186
      PRINT 101,K
  101 FORMAT (/2X,10('*'),2X,'NUMBER OF HISTOGRAMM=',I4,2X,10('*'))
      GO TO (51,52,53,54,55,56,57,58,59,60),K
   51 PRINT 151
      GO TO 61
   52 PRINT 152
      GO TO 61
   53 PRINT 153
      GO TO 61
   54 PRINT 154
      GO TO 61
   55 PRINT 155
      GO TO 61
   56 PRINT 156
      GO TO 61
   57 PRINT 157
      GO TO 61
   58 PRINT 158
      GO TO 61
   59 PRINT 159
      GO TO 61
   60 PRINT 160  
  151 FORMAT (' DISTRIBUTION OF ANGLES FI PERVICH IN MICRORD')
  152 FORMAT (2X,'DISTRIB OF ANGLES FI i-k  IN MICRO RAD') 
  153 FORMAT (2X,'distribution FI  IN MICRORAD')
  154 FORMAT (2X,'DISTRIB moduly fi   IN MICRO RAD')
  155 FORMAT (2X,' distrib i-j for 4 part ')
  156 FORMAT (2X,' distrib fi for 4 part')
  157 FORMAT (2X,'distrib moduly fi for 4 part ')
  158 FORMAT (2X,'ANGLES BY MONTE-CARLO')
  159 FORMAT (2X,'MODULY UULOV  ANGLES BY MONTE-CARLO')
  160 FORMAT (2X,'i-j ugly PO M-K')   
   61 IHIS=IHS(K)
      IK=IHIS+NLM-5
      LI=IHI(IK+3)
      LLKI=IHIS-1
      L=0
      SR(K)=0.0
      SS(K)=0.0
      IF (AXX(K).LT.1.5) GO TO 1
      SR(K)=XX(K)/AXX(K)
      SX2=XX2(K)
      SXX=XX(K)
      ANX=AXX(K)
      CNX=(SX2-((SXX**2)/ANX))/(ANX-1.0)
      IF (CNX.GT.0.0) GO TO 2
      PRINT 204,CNX,SX2,SXX,ANX
  204 FORMAT (2X,'IN SABR PHIST WE HAV SQRT(-)'/2X,4E12.6)
      GO TO 6
    2 SS(K)=SQRT(CNX)
    6 PRINT 202,GR(K),SA(K),AXX(K),SR(K),SS(K)
  202 FORMAT (2X,'LEFT GRAN AND STP OF HIST=',2F12.6/2X,
     *'ALL NUMB IN HIST MIDL VAL AND STANDARD DEV='/2X,
     *F8.1,2X,F12.6,2X,F12.6)
      DO 4 M=IHIS,IK
      L=L+IHI(M)
      LLKI=LLKI+1
      IF (L.EQ.LI) GO TO 5
    4 CONTINUE
    5 IF (K.LT.5) GO TO 72
      IF (K.GT.7) GO TO 72
      M=IHIS
      DO 132 J=1,80
      A=IHI(M)
      EE(J)=A/4.0
      M=M+1
      JJ=J
      IF (M.GT.LLKI) GO TO 84
  132 CONTINUE
   84 PRINT 244,(EE(M),M=1,JJ) 
   72 PRINT 203,(IHI(M),M=IHIS,LLKI)
      IK=IK+1
  203 FORMAT (2X,10I6)
      L=IK+3
      IF (K.LT.5) GO TO 71
      IF (K.GT.7) GO TO 71
      DO 73 M=IK,L  
   73 IHI(M)=IHI(M)/4
   71 PRINT 205,(IHI(M),M=IK,L)
  205 FORMAT (2X,'OVER LEFT GRAN =',I6,2X,'GREATER HIST=',I6/2X,
     *'EVENTS IN HISIOG=',I6,2X,'ALL EVENTS=',I6)
  102 FORMAT (2X,10('*'),2X,'NUMBER OF HISTOGRAMM=',I4,2X,10('*'))
      PRINT 102,K
      IF (K.LT.2) GO TO 1
      IF (K.GT.2) GO TO 1
      A=AXX(3)/AXX(2)
      IE=IHIS
      DO 44 KK=1,40
      EE(KK)=A*IHI(IE)
   44 IE=IE+1
      PRINT 244,(EE(KK),KK=1,40)
  244 FORMAT (10F6.1)   
    1 CONTINUE
      DO 11 J=1,MM
   11 F(J)=(F(J)-SR(3))/SS(3)
      CALL ELEFUN
      RETURN
      END
C 261
      SUBROUTINE ELEFUN
C F- EMPIRIKAL FUNKTION DISTR N- NAMBER OF OBS
C 31 DEK 1998 LEPEKHIN
      COMMON /ELEF/F(2000),N
      DIMENSION P(4)
      DATA P/1.62,1.747,1.862,2.001/
      DATA AMX,AMI/1.0E12,-1.0E12/
      DATA PSKM,PAKK/0.743,1.63/
      CALL RANGIR
      AMAX=AMI
      AMIN=AMX
      AE=N
      C=0.0
      AE2=AE*2.0
      SKM=0.0
      AKK=AMI
      DO 3 K=1,N
      C=C+1.0
      FF=F(K)
      CALL NDTR(FF,PP,DD)
      FF=PP
      BBE=FF-C/AE
      SKM=SKM+(FF-(2.0*C-1.0)/AE2)**2
      IF (BBE.GT.AKK) AKK=BBE
      IF (AMAX.LT.BBE) AMAX=BBE
      IF (AMIN.GT.BBE) AMIN=BBE
    3 CONTINUE
      X=AMAX-AMIN+1.0/AE
      V=X*(SQRT(AE)+0.155+0.24*SQRT(1.0/AE))
      AKK=AKK*SQRT(AE)
      SKM=SKM+1.0/(12.0*AE)
      PRINT 104,SKM,PSKM,AKK,PAKK,V,P(4)
  104 FORMAT (
     *2X,'CRITR KRAM-MIS AND FOR 1 PERCENT',2E10.3/
     *2X,'CRITER KOLMOGOROV AND FOR 1%',2E10.3/
     *2X,'CUIPR CRITER AND FRO 1%',2E10.3)
      RETURN
      END
C 299
      SUBROUTINE NDTR(X,P,D)
      AX=ABS(X)
      T=1.0/(1.0+0.2316419*AX)
      D=0.3989423*EXP(-(X*X)/2.0)
      P=1.0-D*T*((((1.330274*T-1.821256)*T+1.781478)*T-
     *0.3565638)*T+0.3193815)
      IF (X) 1,2,2
    1 P=1.0-P
    2 RETURN
      END
C 311
C 
      SUBROUTINE GAUSS(S,X)
      A=0.0
      DO 1 K=1,12
      CALL CRND(R)
    1 A=A+R
      X=(A-6.0)*S
      RETURN
      END
C
C
      SUBROUTINE CRND(R)
C GENERATION OF RANDOM DISTRIB  NUMBER FROM 1.0 TO 0.0
C 1 DECEMBER LEPEKHIN
      DIMENSION C(50)
      DATA INI,E/0,1.0/
      DATA C/
     *  2.0,3.0,5.0,7.0,11.0,13.0,17.0,19.0,23.0,29.0,31.0,
     *37.0,41.0,43.0,47.0,53.0,59.0,61.0,63.0,71.0,73.0,79.0,83.0,
     *89.0,97.0,101.0,103.0,107.0,109.0,113.0,127.0,131.0,137.0,
     *139.0,149.0,151.0,157.0,163.0,167.0,173.0,179.0,181.0,
     *191.0,193.0,197.0,199.0,211.0,223.0,223.0,229.0/
      IF (INI.GT.1) GO TO 2
      DO 1 K=1,50
    1 C(K)=SQRT(C(K))
      INI=2
    2 I=50.0*RND()+1
      IF (I.GT.50) I=50
      R=RND()+C(I)
      C(I)=R
      R=AMOD(R,E)
      RETURN
      END
C 
      FUNCTION RND()
      DATA I,M,L,A/272667285,2147483646,65539,0.4656613E-9/
    1 J=I*L
      IF (I) 5,6,6
    5 J=J+M+1
    6 I=J
      R=I
      R=R*A
      IF (R.LT.0.0) GO TO 1
      IF (R.GT.1.0) GO TO 1
      RND=R
      RETURN
      END