C	-+-+-+-+-+   H A E S C M . F T   -+-+-+-+-+
C
C	COMPUTE SOME VALUES IN THE X Y AREA
C
C	UP  2-MAR-84  ADDS ANGULAR DISTRBUTION HISTOGRAM
C	UP  4-MAR-84  PLOTS THE ANGULAR DISTRIBUTION HISTOGRAM
C
	SUBROUTINE SCMPUT
	INCLUDE HAEBUF.FI
	INCLUDE HAETTY.FI
	INCLUDE HAEGSA.FI
	INCLUDE HAEPTI.FI
	INCLUDE HAPPEN.FI
	INCLUDE HAECSZ.FI
	INCLUDE HABRK.FI
	REAL XYSAM,TIME,ATAN2,COND
	INTEGER TOASCI,NCOND
	EXTERNAL XYSAM,TIME,TOASCI,NCOND,ATAN2,COND
C
	REAL DELX,DELY,XY
	REAL XAMED,YAMED,		@ ARITHMETISCHER MITTELWERT ALLER PUNKTE
     *  CPUTIM,				@ HERE WE STORE THE COMPUTERS TIME NEEDED
     *  MERADI,DMERAD,			@ MITTLERE GEOMETRISCHE ENTFERNUNG (XAMED,YAMED) ZU ALLEN PUNKTEN 
C					@ DAS IST EIN KREIS UM DEN MITTELPUNKT MIT RADIUS MERADI, 
C					@ DMERAD IST DIE STANDARTABWEICHUNG DAZU
     *  MEXY,				@ MITTLERER WEG ZWISCHEN ZWEI MESSPUNKTEN
     *  DEMEXY,				@ STANDARTABWEICHUNG DAZU
     *  PERIOD,				@ 1/SAMRAT
     *  FLACH				@ MERADI**2 * 3.1416 :  DIE FLAECHE DES KREISES UM XAMED,YAMED
     *  ,ANGHIS(36),WINK(36)		@ ANGULAR DISTRIBUTION HISTOGRAM
	INTEGER I,J,POS1,POSIT,
     *  S0,S1,				@ BEGIN, END OF EXAMINATION
     *  MM,SS,MMM,SSS,TE,		@ SOME TEMPS
     *  CA,CO				@ CHANNEL OF ABZISSE, ORDINATE
C
C	STATEMENT FUNCTIONS:
C
	XSAM(IP1)=XYSAM(IP1,CA)
	YSAM(IP1)=XYSAM(IP1,CO)
	ANGLE(IP3)=MOD(INT((ABS(ATAN2(YSAM(IP3)-YAMED,XSAM(IP3)-XAMED))
     *   +COND(YSAM(IP3)-YAMED.LT.0,3.1416))*5.7297+.5),36)+1 @ THIS REALLY WORKS
	PHYTAG(P1,P2)=SQRT(P1*P1+P2*P2)
	DIST (IP3)   =PHYTAG(XAMED-XSAM(IP3),YAMED-YSAM(IP3))
	XCON(IP4)=COND(ABS(XSAM(IP4)-XSAM(IP4-1)).GT.3,
     *   XSAM(IP4)-XSAM(IP4-1))
	YCON(IP4)=COND(ABS(YSAM(IP4)-YSAM(IP4-1)).GT.3,
     *   YSAM(IP4)-YSAM(IP4-1))
C
	IF (SAMCNT.EQ.0) RETURN		@ NO SAMPLE TOKEN, RETURN
	IF (OPTION.NE.12) CALL ASKHIM(4) @ IF SC* THEN DO NOT ASK FOR BEGIN,SPAN
	IF (BREAK(11)) GOTO 70		@ USER GETS RID OF PROGRAM
	CPUTIM=TIME(I)			@ GET THE TIME OF START
C
	WRITE (3,2) 
	CA=0
	CO=1
	DO 50 POS1=BEGIN*SAMRAT,ENDS*SAMRAT-1,OVRLAP*SAMRAT
	POSIT=POS1+SPAN*SAMRAT		@ THIS IS THE END OF THE ACTUAL SCAN
	IF (POSIT.GT.ENDS*SAMRAT) GOTO 60 @ SKIP LAST OVERLAPPING COMPUTATION
	S0=POS1
	S1=S0+SPAN*SAMRAT
	XAMED=XSAM(S0)
	YAMED=YSAM(S0)
C
	DO 10 I=S0+1,S1,1
	XAMED=XSAM(I)+XAMED
	YAMED=YSAM(I)+YAMED
10	CONTINUE
	YAMED=YAMED/(S1-S0)		@ ARITHMETISCHER MITTELWERT
	XAMED=XAMED/(S1-S0)
C
C	BERECHNE DIE WEGLAENGE
C
	XY=0
	IF(BREAK(11)) GOTO 70
	DO 15 I=S0+1,S1			@ S1-S0 PUNKTE

15	XY=XY+PHYTAG(XCON(I),YCON(I))
C
C	MITTLERE WEGLAENGE IM EINEM MESSINTERVALL (1/SAMRAT SEC.)
C
	MEXY=XY/(S1-S0)			@ MITTLERE WEGLAENGE ZWISCHEN ZWEI PUNKTEN --> MEXY
	DEMEXY=0			@ UND JETZT BERECHNEN WIR DIE STANDARTABWEICHUNG
	DO 17 I=S0+1,S1			@ S1-S0 PUNKTE
17	DEMEXY=(MEXY-PHYTAG(XCON(I),YCON(I)))**2+DEMEXY @ S1-S0-1 PUNKTE NUR!
	DEMEXY=SQRT(DEMEXY/(S1-S0-1))	@ STANDARTABWEICHUNG --> DEMEXY
C
C	MITTLERER KREIS UM (XAMED,YAMED)
C
	MERADI=0
	DO 20 I=S0,S1			@ S1-S0+1 PUNKTE
20	MERADI=MERADI+DIST(I)
	MERADI=MERADI/(1+S1-S0)		@ MITTLERER KREIS
C
C	STANDARTABWEICHUNG ZU MERADI
C
	DMERAD=0
	DO 30 I=S0,S1			@ S1-S0+1 PKT
30	DMERAD=(MERADI-DIST(I))**2+DMERAD
	DMERAD=SQRT(DMERAD/(S1-S0))
C
C	ANGULAR DISTRIBUTION HISTOGRAM ( HUFSCHMIDT ET AL.)
	IF (BREAK(11)) GOTO 70
	CALL MOVE (-36,0,ANGHIS)
	DO 35 I=S0,S1			@ ANOTHER SCAN 
	TE=ANGLE(I)
35	ANGHIS(TE)=ANGHIS(TE)+1		@ HERE WE COMPUTE THE HISTOGRAM
C
C
C	PREPARE THE OUTPUT LINE
C
	MM=TOASCI(S0/SAMRAT/60)
	SS=TOASCI(MOD(S0/SAMRAT,60))
	MMM=TOASCI((S1/SAMRAT)/60)
	SSS=TOASCI(MOD(S1/SAMRAT,60))
C
	PERIOD=1./SAMRAT
	FLACH =MERADI**2*3.141593
	WRITE (3,1) LABEL,MM,SS,MMM,SSS,
     *  XAMED,YAMED,MERADI,FLACH,DMERAD,XY,PERIOD,MEXY,DEMEXY
	WRITE (3,4) ((ANGHIS(I+J),J=1,28,9),I=0,8)
C
C	HERE WE PLOT THE SWAY DIRECTION HISTOGRAM
	DO 110 I=1,36
110	WINK(I)=(I*10.-5.)*.017453	@ FILL 0 TO 2*PI IN 5 DEGREE STEPS
	IF (BREAK(11)) GOTO 70
	CALL TWAIT (10000)		@ WAIT FOR THE PRINTER TO COMPLETE
	CALL STPLT
	CALL XYPLOT (XOFSET,YOFSET,-PENUP)
	CALL XYPLOT (XLEN/2.,(YLEN-1.)/2.,-PENUP) @ SET ORIGIN FOR POLAR PLOT RTN
	CALL POLAR (ANGHIS,WINK,36,1,0,3,0,(YLEN-1.)/200.)
	CALL XYPLOT(-XLEN/2.,-(YLEN-1.)/2.,-PENUP)
	CALL LABPLT(0,-YOFSET,ENDS)
	IF (BREAK(11)) GOTO 120
	CALL SYMBOL (-YZ,2*YZ,YZ,'SWAY DIRECTION HISTOGRAM',90.,25)
	CALL SYMBOL(YZ,YLEN+YZ,YZ,LABEL,0,MIN0(42,LABCNT*6)) @ PLOT 42 CHARACTERS INTO THE FIRST LINE OF LABEL
	IF (LABCNT.GT.7)
     *  CALL SYMBOL (YZ,YLEN,YZ,LABEL(8),0,LABCNT*6-42) @ AND THE REMAINING INTO THE NEXT LINE
120	CALL XYPLOT (38.,25.,PENUP)
	CALL EXPLT
	IF (BREAK(11)) GOTO 70		@ USER GET'S RID OF PLOTTER
C
C
50	CONTINUE
60	CONTINUE
	CPUTIM=TIME(I)-CPUTIM		@ COMPUTE EXECUTION TIME
	WRITE (TTO,3) CPUTIM
70	RETURN
1	FORMAT (1X,10A6,T100,A2,1H:,A2,' --> ',A2,1H:,A2,
     *  /'  MITTE X',F7.1,'  Y',F7.1,5X,
     *  'RADIUS',F6.1,' FLAECHE',F8.0,' S.A.',F5.1,
     *  3X,'WEG',F7.0,3X,'FUER',F6.3,' SEC.  WEG',F6.3,' S.A.',F6.2)
2	FORMAT (/50(2H -)/1X,13A6)
3	FORMAT (' RECHENZEIT:',F6.1,' SEC')
4	FORMAT ('0 HISTOGRAMM DER WINKELVERTEILUNG:  ARCTAN (YI/XI)',
     *   '  POSITION HISTOGRAM',//,22X,
     *   ,' 1. QUADRANT   2. QUADRANT   3. QUADRANT   4. QUADRANT '//
     *  ,2(5(21X,4F12.0/)/))
	END
                                                                                     