**==aa0001.f  
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*********************************
C**       Program SPOTIN        **
C*********************************
C
      INCLUDE 'implicit_none.finc'
C
C
C PURPOSE : program to display a film image for picking out nodal
C           spot positions for autoindexing. Derived from `DISPLAUE'
C           routine of GENLAUE.
C
C     calculate a radial background curve
C     and use a threshold above this for display
C
C AUTHOR : IJC, PAM
C
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
      INCLUDE 'backgr.finc'
      INCLUDE 'sizes.finc'
C
      COMMON /PLOTS/SCAN_START_X,SCAN_START_Y,IDX,IDY
      INTEGER IDX,IDY,SCAN_START_X,SCAN_START_Y
      COMMON /CMBACK/IBGY,IBGV
      INTEGER IBGY(30),IBGV(30)
C
      INCLUDE 'fiducials.finc'
      INCLUDE 'observed_spots.finc'
      INCLUDE 'flmconsts.finc'
C
      COMMON /OMEGA/COSOM,SINOM
      REAL COSOM,SINOM
      COMMON/IPLATE/ IPLATE
      LOGICAL IPLATE
C
C====== FUNCTION TYPES
C
      INTEGER LUNSTI, LUNSTO, LENSTR
      LOGICAL CCPEXS
C
C====== LOCALS:
C
      CHARACTER*20 TRM
      CHARACTER*80 FILM_NAME
      CHARACTER*60 PATH
      CHARACTER*40 NAME
      CHARACTER*12 TYPE
      CHARACTER*6 VERS
      CHARACTER*60 LOGNAM
      CHARACTER*80 FILTMP
      CHARACTER*1 KCH
      INTEGER L, ITTYP
      INTEGER FILM_UNIT, NUM_CURS
      INTEGER IBB(MAX_RASTERS)
      CHARACTER*1 ANS
      CHARACTER*60 IMG_FILE,  FNAME,  GEN_FILE
      LOGICAL QUIT, IMSAV, IMREP, EXISTS, INDSPOT, BEAMC
      LOGICAL DRAW_VEC
      INTEGER I, J, JJ, K, KK, NUM_PLOT_X, NUM_PLOT_Y, NREC, NYREC
      INTEGER NEXT_BACK_SEG, COFG
      INTEGER JTHRESH, ITHRESH
      REAL OMEGA, OMRAD, X, Y
      INTEGER IRECL, IER, IBLK, CUR_BACK_SEG, VEC_BEGIN
      INTEGER IADDR, IADR, IVAL, VEC_END
      INTEGER II, IBACK, IX, IY
      REAL XX, YY
      REAL XM, YM, XP, YP
      INTEGER IRAD
      INTEGER NUM_OVER
      INTEGER IFAIL
      LOGICAL ERR
      DATA FILM_UNIT/3/
C
C-----------------------------------------------------------------------
C
      NUM_PLOT_X = 700
      NUM_PLOT_Y = 700
      MAX_OBSERVED = MAX_OBS
      IN_TERM = LUNSTI(0)
      OUT_TERM = LUNSTO(0)
      IN_DATA = 2
      OUT_DATA = 1
      CALL TRWSET(LUNSTI(0),LUNSTO(0),0,0,0)
      ITTYP = 0
      CALL UGTENV('TERM',TRM)
      IF (TRM(1:5).EQ.'xterm') ITTYP = 1
      CALL QGINIT(0,0,ITTYP,LUNSTI(0),LUNSTO(0))
      CALL TWRITE(' ')
      CALL TWRITE('Program:     SPOTIN')
      CALL TWRITE('Last update: 10/10/91')
      CALL TWRITE('Authors:     I.J.Clifton, P.A.Machin')
      CALL TWRITE(' ')
 20   CALL TPROMP('Enter name of .gen file (default ext=.gen): ')
      CALL TREAD(FNAME,QUIT)
      CALL FILEXT(FNAME,'gen',GEN_FILE,I)
      EXISTS = CCPEXS(GEN_FILE)
      IF (.NOT.EXISTS) GO TO 20
      IFAIL = 1
      CALL CCPDPN(-9,GEN_FILE,'READONLY','F',0,IFAIL)
      IF (IFAIL.EQ.-1) GO TO 20
C
C====== Read in `nearly-constants' from .GEN file or LAUEFLM.DAT
C
      CALL READCONSTS(9)
      CLOSE (9)
C
C====== See if film or image plate
C
      IPLATE = .FALSE.
      CALL TPROMP('Processing film (f) or image plate data (i) [f]: ')
      CALL TREAD(KCH,QUIT)
      IF (KCH.EQ.'i' .OR. KCH.EQ.'I') IPLATE = .TRUE.
C
C====== Initialise raster array
C
      DO 30 I = 1 , MAX_RASTERS
         IBB(I) = 0
 30   CONTINUE
C
C====== Initialise plot variables
C
      IDX = X_RASTERS/NUM_PLOT_X
      IDY = Y_RASTERS/NUM_PLOT_Y
      IF (IDX.NE.IDY) THEN
         NUM_PLOT_X = MIN(X_RASTERS,Y_RASTERS) / MAX(IDX,IDY)
         NUM_PLOT_Y = NUM_PLOT_X
         IDX = MAX(IDX,IDY)
         IDY = IDX
      END IF
      NREC = NUM_PLOT_X*IDX
      NYREC = NUM_PLOT_Y*IDY
      JTHRESH = 50
      IF (IPLATE) JTHRESH = 500
C
C====== Get film filename (assuming defaults)
C
 40   CALL TPROMP('Film file name (default ext=.od): ')
      CALL TREAD(FILTMP,QUIT)
      CALL FILEXT(FILTMP,'od',FILM_NAME,I)
C
C====== Assume default device is film: or from env variable FILM
C
      CALL CCPPSF(FILM_NAME,PATH,NAME,TYPE,VERS)
      IF (PATH.EQ.' ') THEN
         CALL UGTENV('FILM',LOGNAM)
         IF (LOGNAM.NE.' ') THEN
            FILTMP = FILM_NAME
            FILM_NAME = LOGNAM
            L = LENSTR(FILM_NAME)
            FILM_NAME(L+1:) = FILTMP
         END IF
      END IF
C
      EXISTS = CCPEXS(FILM_NAME)
      IF (.NOT.EXISTS) THEN
         CALL TWRITE('!!! File not found')
         GO TO 40
      END IF
C
C====== Open film file
C
      IF (IPLATE) THEN
         IRECL = 2*Y_RASTERS
         CALL IMPOPN(-FILM_UNIT,FILM_NAME,IRECL,SCAN_START_X,OUT_TERM,
     +               ERR)
      ELSE
         IRECL = Y_RASTERS
         CALL FLMOPN(-FILM_UNIT,FILM_NAME,IRECL,SCAN_START_X,OUT_TERM,
     +               ERR)
      END IF
C
C====== Find fiducials
C
      IF (.NOT.IPLATE.AND.NUM_FID.GT.0) CALL FIND_FIDUCIALS(FILM_UNIT)
C
C====== Calculate background
C
      CALL CALC_BACKGR(FILM_UNIT)
      OMEGA = W_F + W_C
      OMRAD = OMEGA*3.141592653/180.
      COSOM = COS(OMRAD)
      SINOM = SIN(OMRAD)
C
      SCAN_START_X = X_CEN_F - NREC/2 + 1
      SCAN_START_Y = Y_CEN_F - NYREC/2 + 1
C
 50   CALL TWRITE('Do you want to select spots by:-')
      CALL TWRITE('1 = Cursor coordinates or')
      CALL TWRITE('2 = Centre of gravity of box centred at cursor')
      COFG = 2
      CALL TASKI('   ',COFG,QUIT)
      IF (COFG.LT.1 .OR. COFG.GT.2) GO TO 50
      BEAMC = .TRUE.
      CALL TASKQ('Cursor in the beam centre',BEAMC,QUIT)
      INDSPOT = .FALSE.
      CALL TASKQ('Supply spot indices',INDSPOT,QUIT)
C
C====== Output message re cursor input
C
      IF (BEAMC) THEN
         CALL TWRITE('Note when plot of film and nodal spots is done,')
         CALL TWRITE(
     +          'use the cursor to FIRST input the BEAM CENTRE and then'
     +          )
         CALL TWRITE('input nodal spots on film. ')
         CALL TWRITE('Type e to finish. ')
      ELSE
         CALL TWRITE('Note when plot of film and nodal spots is done,')
         CALL TWRITE('use the cursor to input nodal spots on film. ')
         CALL TWRITE('Type e to finish. ')
      END IF
C
      IMREP = .TRUE.
      CALL TASKQ('Do you want to replot an image',IMREP,QUIT)
      IF (IMREP) THEN
         CALL TPROMP('Image file name (default ext=.img): ')
         CALL TREAD(FILTMP,QUIT)
         CALL FILEXT(FILTMP,'img',IMG_FILE,I)
         CALL QGENTR
         CALL QGCLR
         CALL QGPIXC
         IF (BEAMC) THEN
            CALL TWRITE('Beam centre first: Type e to finish')
         ELSE
            CALL TWRITE('Type e to finish')
         END IF
         CALL REPLOT(60,IMG_FILE,LUNSTO(0))
         GO TO 110
      END IF
C
      CALL TASKI('Enter above background threshold ',JTHRESH,QUIT)
C
      IMSAV = .FALSE.
      CALL TASKQ('Do you want to save film image',IMSAV,QUIT)
      IF (IMSAV) THEN
         CALL TPROMP('Image file name (default ext=.img): ')
         CALL TREAD(FILTMP,QUIT)
         CALL FILEXT(FILTMP,'img',IMG_FILE,I)
         CALL STASAV(61,IMG_FILE)
      END IF
C
C====== Start of image
C
      IF (IPLATE) THEN
         CALL IMPPOS(FILM_UNIT,SCAN_START_X,OUT_TERM,ERR)
      ELSE
         CALL FLMPOS(FILM_UNIT,SCAN_START_X,OUT_TERM,ERR)
      END IF
C
C====== Initialise plotter package
C
      CALL QGENTR
      CALL QGCLR
      CALL QGPIXC
      IF (BEAMC) THEN
         CALL TWRITE('Beam centre first: Type e to finish')
      ELSE
         CALL TWRITE('Type e to finish')
      END IF
C
C====== Do the display
C
      IBLK = SCAN_START_X
      DO 100 J = 1 , NUM_PLOT_X
C
C====== Set background divisions in y for this scan
C
         CALL BKSET(IBLK)
         DO 90 JJ = 1 , IDX
            ITHRESH = IBGV(1) + JTHRESH
            NEXT_BACK_SEG = 2
            CUR_BACK_SEG = IBGY(1)
            IF (IPLATE) THEN
               CALL IMPLRD(FILM_UNIT,IBB,Y_RASTERS,6,ERR)
            ELSE
               CALL FILMRD(FILM_UNIT,IBB,Y_RASTERS,6,ERR)
            END IF
            IBLK = IBLK + 1
C
C====== Loop over record and plot
C
            DRAW_VEC = .FALSE.
            IADDR = SCAN_START_Y
C
C====== Find start of background
C
 60         IF (IADDR.GE.CUR_BACK_SEG) THEN
               CUR_BACK_SEG = IBGY(NEXT_BACK_SEG)
               ITHRESH = IBGV(NEXT_BACK_SEG) + JTHRESH
               NEXT_BACK_SEG = NEXT_BACK_SEG + 1
               GO TO 60
 
            END IF
 
            DO 80 K = 1 , NUM_PLOT_Y
               IVAL = 0
C
C====== Set ival to max. of the (idy) scan pts at this plot pt.
C
               DO 70 KK = 1 , IDY
                  IF (IBB(IADDR).GT.IVAL) IVAL = IBB(IADDR)
                  IADDR = IADDR + 1
 70            CONTINUE
C
               IF (IADDR.GE.CUR_BACK_SEG) THEN
                  ITHRESH = IBGV(NEXT_BACK_SEG) + JTHRESH
                  CUR_BACK_SEG = IBGY(NEXT_BACK_SEG)
                  NEXT_BACK_SEG = NEXT_BACK_SEG + 1
               END IF
C
               IF (IVAL.GT.ITHRESH .AND. .NOT.DRAW_VEC) THEN
                  VEC_BEGIN = K
                  DRAW_VEC = .TRUE.
                  CALL QGMOVE(FLOAT(J),FLOAT(VEC_BEGIN))
                  CALL PUTSAV
               END IF
C
               IF (IVAL.LE.ITHRESH .AND. DRAW_VEC) THEN
                  VEC_END = (K-1)
                  DRAW_VEC = .FALSE.
                  CALL QGDRAW(FLOAT(J),FLOAT(VEC_END))
                  CALL PUTSAV
               END IF
 
 80         CONTINUE
C
            IF (DRAW_VEC) THEN
               CALL QGDRAW(FLOAT(J),FLOAT(NUM_PLOT_Y))
               CALL PUTSAV
            END IF
 90      CONTINUE
C
 100  CONTINUE
C
C====== End of film plot
C
      IF (IMSAV) CALL FINSAV
C
C====== Cursor in [beam centre and] spots
C
 110  IF (BEAMC) THEN
         CALL QGCURS(XM,YM,KCH)
         CALL QGMPXY(XM,YM,XP,YP)
         IX = NINT(XP)
         IY = NINT(YP)
         IF (COFG.EQ.2) THEN
C
C====== Box input mode.
C
C====== Calculate background using radial curve and spot position
C====== change to scanner units
C
            XX = IX*IDX + SCAN_START_X
            YY = IY*IDY + SCAN_START_Y
C
C====== NB Background is centred on (xcen,ycen)
C
            IRAD = SQRT((XX-XCEN)**2+(YY-YCEN)**2)
            II = N_BACKGR - (IRAD/N_PER_BIN)
 
C
C====== Overestimate background
C
            IBACK = BACKGR(II) + 10
            IX = XX
            IY = YY
            CALL C_OF_G(FILM_UNIT,IX,IY,20,20,IBACK,X,Y,NUM_OVER)
            IF (NUM_OVER.LE.0) THEN
               CALL TWRITE('!!! no pixels above background level')
               STOP
 
            END IF
C
         ELSE
C
C====== Simple input mode.
C
C====== Change to scanner units
C
            X = IX*IDX + SCAN_START_X
            Y = IY*IDY + SCAN_START_Y
         END IF
         X_C = X - X_CEN_F
         Y_C = Y - Y_CEN_F
      END IF
      CALL NODEIN(NUM_CURS,COFG,INDSPOT,FILM_UNIT)
C
      IF (IPLATE) THEN
         CALL IMPCLS(FILM_UNIT)
      ELSE
         CALL FLMCLS(FILM_UNIT)
      END IF
      CALL QGCLOS
      END
**==bkset.f  
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*************************
C**       BKSET         **
C*************************
C
C
      SUBROUTINE BKSET(IX)
C
      INCLUDE 'implicit_none.finc'
C
C Purpose : For a given record set yvalues where background
C           changes (in ibgy) and store the background to
C           use in ibgv. note there are usually an odd
C           number of regions, except when ix=x_cen_f
C
      INTEGER IX
C
      INCLUDE 'backgr.finc'
C
      COMMON /CMBACK/IBGY,IBGV
      INTEGER IBGY(30),IBGV(30)
C
      INCLUDE 'fiducials.finc'
C
      COMMON /PLOTS/SCAN_START_X,SCAN_START_Y,IDX,IDY
      INTEGER IDX,IDY,SCAN_START_X,SCAN_START_Y
C
C Author: IJC
C
C====== LOCALS:
C
      INTEGER RAD,RAD_MAX,Y_BIN,I,NUM_SEEN
      REAL Y_DIFF
C
C-------------------------------------------------------------------------------
C
      RAD_MAX = Y_RANGE/2
      DO 20 I = 1 , 30
         IBGY(I) = 10000
         IBGV(I) = BACKGR(1)
 20   CONTINUE
      RAD = RAD_MAX - N_PER_BIN
      I = 0
C
C====== See if we cut the next background region
C
 30   IF (IABS(IX-XCEN).LT.RAD) THEN
         I = I + 1
         Y_DIFF = RAD**2 - (IX-XCEN)**2
         Y_BIN = YCEN - SQRT(Y_DIFF)
         IBGY(I) = Y_BIN
         IBGV(I) = BACKGR(I)
C
C====== Seen all the regions?
C
         IF (I.LT.N_BACKGR) THEN
            RAD = RAD - N_PER_BIN
            GO TO 30
         END IF
      END IF
 
      NUM_SEEN = I
      IF (NUM_SEEN.LT.1) RETURN
C
      IBGV(NUM_SEEN+1) = BACKGR(I+1)
C
C====== Apply symmetry
C
      DO 40 I = 1 , NUM_SEEN
         IBGV(NUM_SEEN+I+1) = IBGV(NUM_SEEN-I+1)
         IBGY(NUM_SEEN+I) = 2*YCEN - IBGY(NUM_SEEN-I+1)
 40   CONTINUE
      RETURN
      END
**==c_of_g.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*************************
C**       C_OF_G        **
C*************************
C
C
      SUBROUTINE  C_OF_G (FILM_UNIT, XC, YC, X_SIZE, Y_SIZE, THRESH,
     +                    CGX, CGY, NUM_OVER)
C
C Purpose: Finds c of g of intensities lying in box centred on (xc, yc)
C          of dimensions 2*x_size+1 and 2*y_size+1,
C          allowing for a constant intensity threshold
C
      INCLUDE 'implicit_none.finc'
C
      INTEGER    FILM_UNIT
C
C====== Integers because in SCANNER space
C
      INTEGER    XC, YC, X_SIZE, Y_SIZE
      INTEGER    THRESH, NUM_OVER
      REAL       CGX, CGY
C
C====== NUM_OVER: number of spots in the box over the threshold.
C
      INCLUDE 'sizes.finc'
C
      COMMON/IPLATE/ IPLATE
      LOGICAL IPLATE
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
C====== LOCALS:
C
      LOGICAL    ERR
      CHARACTER*80 STR
      INTEGER    IX, IY, NSTRIP
      INTEGER    ISTRIP(MAX_RASTERS)
      INTEGER    INTENS, I
      REAL       SUMX, SUMY, SUM_INTENS, DX, DY
C
C------------------------------------------------------------------------------
C
      NUM_OVER = 0
      SUMX = 0.0
      SUMY = 0.0
      SUM_INTENS = 0.0
C
      NSTRIP = XC - X_SIZE
      IF (XC-X_SIZE.LT.1 .OR. XC+X_SIZE.GT.X_RASTERS .OR. 
     +    YC-Y_SIZE.LT.1 .OR. YC+Y_SIZE.GT.Y_RASTERS) THEN
         WRITE (STR,6010) XC , YC
         CALL TWRITE(STR)
         NUM_OVER = 0
         CGX = 0.0
         CGY = 0.0
         RETURN
      END IF
      IF (IPLATE) THEN
         CALL IMPPOS(FILM_UNIT,NSTRIP,OUT_TERM,ERR)
      ELSE
         CALL FLMPOS(FILM_UNIT,NSTRIP,OUT_TERM,ERR)
      END IF
C
      DO 30 IX = -X_SIZE , X_SIZE
         IF (IPLATE) THEN
            CALL IMPLRD(FILM_UNIT,ISTRIP,Y_RASTERS,OUT_TERM,ERR)
         ELSE
            CALL FILMRD(FILM_UNIT,ISTRIP,Y_RASTERS,OUT_TERM,ERR)
         END IF
         NSTRIP = NSTRIP + 1
C
         DO 20 IY = -Y_SIZE , Y_SIZE
            INTENS = ISTRIP(YC+IY)
            INTENS = MAX(0,INTENS-THRESH)
            IF (INTENS.GT.0) NUM_OVER = NUM_OVER + 1
            SUMX = SUMX + IX*INTENS
            SUMY = SUMY + IY*INTENS
            SUM_INTENS = SUM_INTENS + INTENS
 20      CONTINUE
 
 30   CONTINUE
C
      IF (SUM_INTENS.EQ.0.0) THEN
         DX = 0
         DY = 0
      ELSE
         DX = SUMX/SUM_INTENS
         DY = SUMY/SUM_INTENS
      END IF
C
      CGX = XC + DX
      CGY = YC + DY
C
      RETURN
C
C====== Format statements
C
 6010 FORMAT ('**Attempt to read beyond scanned film** xc =',I5,' yc =',
     +        I5)
      END
**==calc_backgr.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C******************************
C**       CALC_BACKGR        **
C******************************
C
C
      SUBROUTINE  CALC_BACKGR (FILM_UNIT)
C
C Purpose: Estimates a radial background curve
C
      INCLUDE 'implicit_none.finc'
C
      INTEGER    FILM_UNIT
C
      INCLUDE 'backgr.finc'
      INCLUDE 'sizes.finc'
C
      COMMON/IPLATE/ IPLATE
      LOGICAL IPLATE
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
C====== LOCALS:
C
      LOGICAL    ERR
      CHARACTER*80 STR
      INTEGER    I, FREQ(1000), IVAL(1000), NFREQ, IBIN
      INTEGER    X, Y
      INTEGER    BIN, MAX_IN_BIN
      INTEGER    N_INTENS, INTENS, SUM_INTENS(N_BINS), SUM
      INTEGER NUM_PLOT_X, NUM_PLOT_Y, IDY
      PARAMETER (NUM_PLOT_X=700,NUM_PLOT_Y=700)
      INTEGER    ISTRIP(MAX_RASTERS)
C
C------------------------------------------------------------------------------
C
C====== Calculate the approx centre used for background estimation
C
      XCEN = X_RASTERS/2
      YCEN = Y_RASTERS/2
      NFREQ = 256
      IF (IPLATE) NFREQ = 1000
C
      IDY = Y_RASTERS/NUM_PLOT_Y
      Y_RANGE = IDY*NUM_PLOT_Y
      Y_START = YCEN - Y_RANGE/2
      N_PER_BIN = Y_RANGE/N_BINS
C
      CALL TWRITE('=== Estimating background...')
C
      MAX_IN_BIN = 0.7*N_PER_BIN + 1
      DO 20 BIN = 1 , N_BINS
         SUM_INTENS(BIN) = 0
 20   CONTINUE
C
      IF (IPLATE) THEN
         CALL IMPPOS(FILM_UNIT,XCEN-20,OUT_TERM,ERR)
      ELSE
         CALL FLMPOS(FILM_UNIT,XCEN-20,OUT_TERM,ERR)
      END IF
 
      DO 80 X = XCEN - 20 , XCEN + 20
         IF (IPLATE) THEN
            CALL IMPLRD(FILM_UNIT,ISTRIP,Y_RASTERS,OUT_TERM,ERR)
         ELSE
            CALL FILMRD(FILM_UNIT,ISTRIP,Y_RASTERS,OUT_TERM,ERR)
         END IF
         Y = Y_START
C
         DO 70 BIN = 1 , N_BINS
            DO 30 I = 1 , NFREQ
               FREQ(I) = 0
               IVAL(I) = 0
 30         CONTINUE
            DO 40 I = 1 , N_PER_BIN
               INTENS = ISTRIP(Y)
               IF (IPLATE) THEN
                  IBIN = INTENS/1000 + 1
                  IF (IBIN.GT.NFREQ) IBIN = NFREQ
               ELSE
                  IBIN = INTENS + 1
               END IF
               FREQ(IBIN) = FREQ(IBIN) + 1
               IVAL(IBIN) = IVAL(IBIN) + INTENS
               Y = Y + 1
 40         CONTINUE
            N_INTENS = 0
            SUM = 0
            DO 50 I = 1 , NFREQ
               N_INTENS = N_INTENS + FREQ(I)
               SUM = SUM + IVAL(I)
               IF (N_INTENS.GT.MAX_IN_BIN) GO TO 60
 50         CONTINUE
 60         SUM_INTENS(BIN) = SUM_INTENS(BIN) + SUM/N_INTENS
 70      CONTINUE
 80   CONTINUE
C
      DO 90 BIN = 1 , N_BINS
C
C====== 41 is range -20..20
C
         SUM_INTENS(BIN) = SUM_INTENS(BIN)/41
 90   CONTINUE
C
C====== Adjust edges
C
      DO 100 I = 1 , N_BINS/2
         BACKGR(I) = MAX(SUM_INTENS(I),SUM_INTENS(N_BINS+1-I))
 100  CONTINUE
 
      BACKGR(1) = BACKGR(2)
      I = N_BINS/2
      BACKGR(I) = BACKGR(I-1) + IABS(BACKGR(I-1)-BACKGR(I-2))
 
      WRITE (STR,6010) (BACKGR(I),I=1,N_BINS/2)
      CALL TWRITE(STR)
      RETURN
C
C====== Format statements
C
 6010 FORMAT (20I4)
      END
**==dump_spots.f  
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*****************************
C**       DUMP_SPOTS        **
C*****************************
C
C
      SUBROUTINE DUMP_SPOTS (NUM_CURS)
C
      INCLUDE 'implicit_none.finc'
C
C PURPOSE : write out the spots cursored in by the user
C NB this version of DUMP_SPOTS applies x_c, y_c and omega
C
      INTEGER NUM_CURS
C
C AUTHOR : IJC
C
      INCLUDE 'observed_spots.finc'
      INCLUDE 'fiducials.finc'
      INCLUDE 'gen_file_header.finc'
      INCLUDE 'sizes.finc'
C
      COMMON /OMEGA/COSOM,SINOM
      REAL COSOM,SINOM
C
C====== LOCALS:
C
      LOGICAL QUIT
      INTEGER I, DUMP_FIL
      REAL X, Y, XX, YY
      PARAMETER (DUMP_FIL=15)
      CHARACTER*60 DUMP_FILE, FILNAM
C
C-------------------------------------------------------------------------------
C
      CALL TPROMP('Dump file name (default ext=.spots): ')
      CALL TREAD(FILNAM,QUIT)
      CALL FILEXT(FILNAM,'spots',DUMP_FILE,I)
      OPEN (FILE=DUMP_FILE,STATUS='NEW',UNIT=DUMP_FIL)
      WRITE (DUMP_FIL,6010) TITLE , CRYSTAL_ID
      WRITE (DUMP_FIL,6020) NUM_CURS
      DO 20 I = 1 , NUM_CURS
C
C====== Convert from scanner units to mm.s
C
         X = (XM_OBS(I)-X_CEN_F-X_C)/MM_RAST_X
         Y = (YM_OBS(I)-Y_CEN_F-Y_C)/MM_RAST_Y
         XX = COSOM*X - SINOM*Y
         YY = COSOM*Y + SINOM*X
         WRITE (DUMP_FIL,6030) NODAL_H_S(I) , NODAL_K_S(I) , 
     +                         NODAL_L_S(I) , XX , YY
 20   CONTINUE
      CLOSE (DUMP_FIL)
      RETURN
C
C======Format statements
C
 6010 FORMAT (1X,A76,/,1X,A6)
 6020 FORMAT (' Number=',I4)
 6030 FORMAT (' h k l  ',3I4,T40,'x y  ',2F6.1)
      END
**==find_fid.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C***************************
C**       FIND_FID        **
C***************************
C
C
      SUBROUTINE  FIND_FID (FILM_UNIT, ROUGH_X, ROUGH_Y, FIDX, FIDY)
C
C Purpose: This routine is supposed to locate a fiducial accurately given
C          its rough position
C
      INCLUDE 'implicit_none.finc'
C
      INTEGER  FILM_UNIT
      REAL     ROUGH_X, ROUGH_Y
      REAL     FIDX, FIDY
C
C Note: if the routine fails to locate the fiducial within 5mm then
C       fidx and fidy are both -1
C
C
      INCLUDE 'sizes.finc'
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
C====== LOCALS:
C
      CHARACTER*80 STR
      REAL     MAX_X_SIZE, MAX_SHIFT, MAX_SHIFT2
C
      INTEGER  XC, YC, X_SIZE, Y_SIZE
      INTEGER  THRESH, NUM_OVER
      REAL     CGX, CGY
      LOGICAL    TRY_AGAIN, QUIT
C
C------------------------------------------------------------------------------
C
C====== Max size is twice starting size
C====== NB fid box is actually fid_box_size/2
C
      MAX_X_SIZE = FID_BOX_SIZE*MM_RAST_X
      MAX_SHIFT = FID_BOX_SIZE/4.*MM_RAST_X
      MAX_SHIFT2 = MAX_SHIFT*MAX_SHIFT
      THRESH = 125
      TRY_AGAIN = .FALSE.
C
      XC = ROUGH_X
      YC = ROUGH_Y
 20   X_SIZE = FID_BOX_SIZE/2.*MM_RAST_X
      Y_SIZE = X_SIZE
C
C====== Assume background is 25 - there is no reasoning behind this value
C
 30   CALL C_OF_G(FILM_UNIT,XC,YC,X_SIZE,Y_SIZE,THRESH,CGX,CGY,NUM_OVER)
C
C====== If flat film, then enlarge box slighty
C====== if box gets to 10mm x 10mm then abort
C====== = 200 x 200
C
      IF (NUM_OVER.EQ.0) THEN
C
         IF (.NOT.TRY_AGAIN) THEN
            TRY_AGAIN = .TRUE.
            CALL TASKQ(
     +   'Couldn''t find fiducial - try with bigger box/lower threshold'
     +   ,TRY_AGAIN,QUIT)
            IF (.NOT.TRY_AGAIN) THEN
               WRITE (STR,6010) ROUGH_X , ROUGH_Y
               CALL TWRITE(STR)
               FIDX = ROUGH_X
               FIDY = ROUGH_Y
               RETURN
            END IF
         END IF
         IF (X_SIZE.LE.MAX_X_SIZE) THEN
C
C====== OK. - expAND the box by 20% - and make sure it does get bigger
C======       else infinite loop.
C
            X_SIZE = X_SIZE*1.2 + 1
            Y_SIZE = X_SIZE
            GO TO 30
C
C====== Can`t find spot in big box - try lowering threshold (weak fiducial?)
C
         ELSE IF (THRESH.GT.75) THEN
            THRESH = THRESH - 25
            WRITE (STR,6020) THRESH
            CALL TWRITE(STR)
            GO TO 30
         ELSE
            WRITE (STR,6030) THRESH
            CALL TWRITE(STR)
            THRESH = THRESH - 25
            CALL TASKI('Enter threshold or 0 to give up',THRESH,QUIT)
            IF (THRESH.GT.0) GO TO 30
C
C====== Failed to locate fiducial - give up!
C
            FIDX = -1
            FIDY = -1
            RETURN
C
         END IF
C
      END IF
C
C====== If we have moved far from the rough place, re-position the box
C
      IF (((CGX-XC)**2+(CGY-YC)**2).GT.MAX_SHIFT2) THEN
C
         XC = CGX
         YC = CGY
         GO TO 20
C
      END IF
C
C====== OK found a spot that hasn't shifted far from box centre
C
      FIDX = CGX
      FIDY = CGY
C
      RETURN
C
C====== Format statements
C
 6010 FORMAT ('Using ESTIMATED fiducial pos. of ',2F8.1)
 6020 FORMAT ('Could not find fiducial - lowering',' threshold to ',I3)
 6030 FORMAT ('STILL no density found above ',I4)
      END
**==find_fiducials.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*********************************
C**       FIND_FIDUCIALS        **
C*********************************
C
C
      SUBROUTINE  FIND_FIDUCIALS (FILM_UNIT)
C
C Purpose: Locates fiducials on film and determines fiducial parameters
C          i.e. x_cen_f, y_cen_f and w_f
C
      INCLUDE 'implicit_none.finc'
C
      INTEGER    FILM_UNIT
C
      INCLUDE 'fiducials.finc'
      INCLUDE 'sizes.finc'
C
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
C====== LOCALS:
C
      CHARACTER*80 STR
      INTEGER    I
      REAL       X_CENTRE, Y_CENTRE
      REAL       RAD_TO_DEG
      PARAMETER (RAD_TO_DEG = 360./(2.*3.1415926))
C
C------------------------------------------------------------------------------
C
      WRITE (STR,6030)
      CALL TWRITE(STR)
C
      DO 20 I = 1 , 3
C
         CALL FIND_FID(FILM_UNIT,FID_POS_X(I),FID_POS_Y(I),FIDX(I),
     +                 FIDY(I))
         IF (FIDX(I).EQ.-1.0) THEN
            WRITE (STR,6100) I
            CALL TWRITE(STR)
            STOP
         END IF
C
 20   CONTINUE
C
      WRITE (STR,6010) FIDX(1) , FIDY(1)
      CALL TWRITE(STR)
      WRITE (STR,6020) FIDX(2) , FIDY(2)
      CALL TWRITE(STR)
      WRITE (STR,6020) FIDX(3) , FIDY(3)
      CALL TWRITE(STR)
C
      IF (X_CEN_F.LE.0) X_CEN_F = (FIDX(1)+FIDX(3))/2.0
      IF (Y_CEN_F.LE.0) Y_CEN_F = (FIDY(1)+FIDY(3))/2.0
C
      W_F = (FIDY(3)-FIDY(2))/(FIDX(3)-FIDX(2))
      W_F = ATAN(W_F)*RAD_TO_DEG
C
      WRITE (STR,6080) RASTER_SIZE*1000. , X_CEN_F , Y_CEN_F
      CALL TWRITE(STR)
      WRITE (STR,6090) W_F
      CALL TWRITE(STR)
C
C====== See if we can find centre by using same procedure for
C====== finding fiducials
C====== Only do this if x_c and y_c are zero
C
      IF (X_C.EQ.0.0 .AND. Y_C.EQ.0.0) THEN
 
         CALL FIND_FID(FILM_UNIT,X_CEN_F,Y_CEN_F,X_CENTRE,Y_CENTRE)
         IF (X_CENTRE.GT.0) THEN
            WRITE (STR,6040) X_CENTRE , Y_CENTRE
            CALL TWRITE(STR)
            X_C = X_CENTRE - X_CEN_F
            Y_C = Y_CENTRE - Y_CEN_F
            WRITE (STR,6050)
            CALL TWRITE(STR)
            WRITE (STR,6060) RASTER_SIZE*1000.
            CALL TWRITE(STR)
            WRITE (STR,6070) X_C , Y_C
            CALL TWRITE(STR)
         ELSE
            WRITE (STR,6110)
            CALL TWRITE(STR)
         END IF
      END IF
      RETURN
C
C====== Format statements
C
 6010 FORMAT ('Fiducials found at:  (',F7.1,',',F7.1,')')
 6020 FORMAT ('                     (',F7.1,',',F7.1,')')
 6030 FORMAT ('Locating fiducials...')
 6040 FORMAT ('Beam centre at:      (',F7.1,',',F7.1,')')
 6050 FORMAT ('Fiducial constants as determined from beam centre:')
 6060 FORMAT ('(units are in ',F5.1,' micron scanner steps)')
 6070 FORMAT ('x_c =',G14.6,',    y_c =',G14.6)
 6080 FORMAT ('Fiducial centre (',F5.1,'u units) ',2G14.6)
 6090 FORMAT ('Omega as determined from fiducials (degrees)',G14.6)
 6100 FORMAT ('!!! Failed to locate fiducial ',I1)
 6110 FORMAT ('!!! Failed to locate beam centre')
      END
**==get_raster_size.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C**********************************
C**       GET_RASTER_SIZE        **
C**********************************
C
C
      SUBROUTINE  GET_RASTER_SIZE(RASTER_SIZE)
C
C Purpose: Match value of raser_size on input against the range of values
C          supported, and then give the user a change to over-ride.
C
      INCLUDE 'implicit_none.finc'
C
      INTEGER NUM_POSS_VALS
      PARAMETER (NUM_POSS_VALS=4)
      REAL RASTER_SIZE
      REAL POSSIBLE_VALS(NUM_POSS_VALS)
      REAL RANGE_LO(NUM_POSS_VALS), RANGE_HI(NUM_POSS_VALS)
      INTEGER I
      LOGICAL QUIT
      COMMON/IPLATE/ IPLATE
      LOGICAL IPLATE
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      DATA POSSIBLE_VALS/12.5E-3,25.E-3,50.E-3,100.E-3/
      DATA RANGE_LO/10.E-3,20.E-3,40.E-3,90.E-3/,
     + RANGE_HI/15.E-3,30.E-3,60.E-3,110.E-3/
C
C------------------------------------------------------------------------------
C
      IF (.NOT.IPLATE) THEN
         DO 20 I = 1 , NUM_POSS_VALS
            IF (RASTER_SIZE.GE.RANGE_LO(I) .AND. 
     +          RASTER_SIZE.LE.RANGE_HI(I)) THEN
               RASTER_SIZE = POSSIBLE_VALS(I)
               GO TO 30
            END IF
 20      CONTINUE
         CALL TWRITE('Raster size calculated from film dimensions')
         CALL TWRITE('doesn''t correspond to a standard size')
         CALL QSBELL
      END IF
C
C======  Convert to um units during interaction
C
 30   RASTER_SIZE = RASTER_SIZE*1000.
      CALL TASKF('Enter raster size (microns)',RASTER_SIZE,1,QUIT)
      RASTER_SIZE = RASTER_SIZE/1000.
      RETURN
      END
**==nodein.f  
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*************************
C**       NODEIN        **
C*************************
C
C
      SUBROUTINE NODEIN(NUM_CURS, COFG, INDSPOT, FILM_UNIT)
C
      INCLUDE 'implicit_none.finc'
C
C
C PURPOSE : use the cursor to input positions
C           of the spots on the film (for x,y)
C           output header info then spot positions for AUTOindexing
C
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER NUM_CURS, COFG
      LOGICAL INDSPOT
C
      INCLUDE 'backgr.finc'
      INCLUDE 'observed_spots.finc'
      INCLUDE 'fiducials.finc'
C
      COMMON /PLOTS/SCAN_START_X, SCAN_START_Y, IDX, IDY
      INTEGER FILM_UNIT, IDX, IDY, SCAN_START_X, SCAN_START_Y
C
C AUTHOR : PAM
C
C
C====== LOCALS
C
      CHARACTER*1 KCH
      LOGICAL QUIT
      INTEGER I, II, IBACK, IX, IY
      REAL XX, YY
      REAL XM, YM, XP, YP
      INTEGER DIF, DELX, DELY, IRAD, H, K, L
      INTEGER NUM_OVER
      DATA DIF/5/
C
C-----------------------------------------------------------------------
C
      NUM_CURS = 1
C
C====== Input observed spot position
C
 20   CALL QGCURS(XM,YM,KCH)
      CALL QGMPXY(XM,YM,XP,YP)
      IX = NINT(XP)
      IY = NINT(YP)
      IF (KCH.EQ.'e' .OR. KCH.EQ.'E') THEN
C
C====== End of input
C
         CALL QGBELL
         CALL TREAD(KCH,QUIT)
         CALL QGALPH
C
C====== Dump of spot positions
C
         CALL DUMP_SPOTS(NUM_CURS-1)
         RETURN
      ELSE
         IF (COFG.EQ.2) THEN
C
C====== Box input mode.
C
C====== Calculate background using radial curve and spot position
C====== change to scanner units
C
            XX = IX*IDX + SCAN_START_X
            YY = IY*IDY + SCAN_START_Y
C
C====== NB Background is centred on (xcen,ycen)
C
            IRAD = SQRT((XX-XCEN)**2+(YY-YCEN)**2)
            II = N_BACKGR - (IRAD/N_PER_BIN)
 
C
C====== Overestimate background
C
            IBACK = BACKGR(II) + 10
            IX = XX
            IY = YY
            CALL C_OF_G(FILM_UNIT,IX,IY,20,20,IBACK,XM_OBS(NUM_CURS),
     +                  YM_OBS(NUM_CURS),NUM_OVER)
            IF (NUM_OVER.LE.0) THEN
               CALL TWRITE('Flat film - try again')
               GO TO 20
 
            END IF
C
         ELSE
C
C====== Simple input mode.
C
C====== Change to scanner units
C
            XM_OBS(NUM_CURS) = IX*IDX + SCAN_START_X
            YM_OBS(NUM_CURS) = IY*IDY + SCAN_START_Y
         END IF
         IF (INDSPOT) THEN
            CALL TPROMP(
     +  '                                                     Indices: '
     +  )
            READ (IN_TERM,*) H , K , L
            NODAL_H_S(NUM_CURS) = H
            NODAL_K_S(NUM_CURS) = K
            NODAL_L_S(NUM_CURS) = L
         END IF
C
C====== Spot OK - inc count of spots
C
 
         NUM_CURS = NUM_CURS + 1
         GO TO 20
      END IF
 
      END
**==readconsts.f  
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*****************************
C**       READCONSTS        **
C*****************************
C
C
      SUBROUTINE  READCONSTS(IUNIT)
C
C Purpose: Reads `nearly-constant' parameters from the .GEN file
C          (on unit IUNIT) or alternatively `LAUEFLM(.DAT)'.
C          Don't try .GEN file if IUNIT<=0.
C          NB These params are also used by INTLAUE and so some aren't
C          needed in SPOTIN.
C
      INCLUDE 'implicit_none.finc'
C
      INCLUDE 'sizes.finc'
      INCLUDE 'fiducials.finc'
      INCLUDE 'flmconsts.finc'
      INCLUDE 'film_params.finc'
      INCLUDE 'gen_file_header.finc'
C
      INTEGER IUNIT
C
C====== LOCALS:
C
      INTEGER IFAIL
      INTEGER I, XCEN, YCEN
      CHARACTER*80 DUMMY, STR
      REAL VALUES(6)
      INTEGER NVALS
      REAL RASTSIZ
      COMMON/IPLATE/ IPLATE
      LOGICAL IPLATE
      COMMON /IO_STREAMS/ IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
      INTEGER             IN_TERM, OUT_TERM, IN_DATA, OUT_DATA
C
C------------------------------------------------------------------------------
C
C====== Just try LAUEFLM.DAT if unit no. <= 0
C
      IF (IUNIT.GT.0) THEN
C
C====== First try and read from .GEN file - skip most data in recs 3-11
C
         READ (IUNIT,6010) TITLE
         READ (IUNIT,6010) CRYSTAL_ID
         DO 20 I = 1 , 7
            READ (IUNIT,6010) DUMMY
 20      CONTINUE
         READ (IUNIT,FMT=*) Y_SCALE
         READ (IUNIT,6010) DUMMY
         READ (IUNIT,END=30,FMT=*) NUM_FID , 
     +                             (FID_POS_X10U(I),FID_POS_Y10U(I),I=1,
     +                             3)
         READ (IUNIT,*) FILM_WIDTH , FID_BOX_SIZE , N1OD , BASEOD , 
     +                  G1OD , CURVE
         READ (IUNIT,'(A)') DUMMY
         CALL STRVAL(DUMMY,' ,',6,VALUES,NVALS)
         X_CEN_F = 0.0
         Y_CEN_F = 0.0
         RASTSIZ = 0.0
         IF (NVALS.GE.1) X_RASTERS = NINT(VALUES(1))
         IF (NVALS.GE.2) Y_RASTERS = NINT(VALUES(2))
         IF (NVALS.GE.3) FILMPLOT_SCALE = VALUES(3)
         IF (NVALS.GE.4) X_CEN_F = VALUES(4)
         IF (NVALS.GE.5) Y_CEN_F = VALUES(5)
         IF (NVALS.GE.6) RASTSIZ = VALUES(6)
         GO TO 40
      END IF
C
C====== Data missing from .GEN file - try to read LAUEFLM(.DAT)
C
 30   WRITE (STR,6020)
      CALL TWRITE(STR)
      IFAIL = 0
      CALL CCPDPN(IN_DATA,'LAUEFLM','READONLY','F',0,IFAIL)
      READ (IN_DATA,*) NUM_FID , (FID_POS_X10U(I),FID_POS_Y10U(I),I=1,3)
      READ (IN_DATA,*) FILM_WIDTH , FID_BOX_SIZE , N1OD , BASEOD , 
     +                 G1OD , CURVE
      READ (IN_DATA,'(A)') DUMMY
      CALL STRVAL(DUMMY,' ,',6,VALUES,NVALS)
      X_CEN_F = 0.0
      Y_CEN_F = 0.0
      RASTSIZ = 0.0
      IF (NVALS.GE.1) X_RASTERS = NINT(VALUES(1))
      IF (NVALS.GE.2) Y_RASTERS = NINT(VALUES(2))
      IF (NVALS.GE.3) FILMPLOT_SCALE = VALUES(3)
      IF (NVALS.GE.4) X_CEN_F = VALUES(4)
      IF (NVALS.GE.5) Y_CEN_F = VALUES(5)
      IF (NVALS.GE.6) RASTSIZ = VALUES(6)
      CLOSE (IN_DATA)
C
 40   IF (NUM_FID.NE.0 .AND. NUM_FID.NE.3) THEN
         CALL TWRITE('!!! Spotin can only handle 3 fiducials')
         STOP
      END IF
      CENT_IN = .FALSE.
      IF (X_CEN_F.GT.0.0 .AND. Y_CEN_F.GT.0.0) CENT_IN = .TRUE.
C
C====== Default for film width is 125mm (standard film)
C
      RASTER_SIZE = RASTSIZ/1000.0
      IF (FILM_WIDTH.EQ.0.0) FILM_WIDTH = 125.0
      IF (RASTER_SIZE.LE.0) THEN
         RASTER_SIZE = FILM_WIDTH/X_RASTERS
C
C====== Set raster size to the nearest value that SPOTIN supports.
C
         CALL GET_RASTER_SIZE(RASTER_SIZE)
         MM_RAST_X = 1./RASTER_SIZE
         MM_RAST_Y = MM_RAST_X*Y_SCALE
      ELSE
         MM_RAST_X = 1./RASTER_SIZE
         MM_RAST_Y = MM_RAST_X*Y_SCALE
      END IF
C
C====== Approx centre of film
C
      XCEN = X_RASTERS/2
      YCEN = Y_RASTERS/2
      IF (NUM_FID.EQ.0) THEN
         XCEN = X_CEN_F
         YCEN = Y_CEN_F
      END IF
C
C====== Now convert fid posns from 10u units to rasters.
C
      DO 50 I = 1 , 3
         FID_POS_X(I) = (FID_POS_X10U(I)/100.*MM_RAST_X) + XCEN
         FID_POS_Y(I) = (FID_POS_Y10U(I)/100.*MM_RAST_Y) + YCEN
 50   CONTINUE
      REWIND (IUNIT)
C
C====== Format statements
C
 6010 FORMAT (A)
 6020 FORMAT ('Warning - data missing in .gen file'/
     +        '  -trying to read from LAUEFLM')
      END
