	PROGRAM TRACE
C
C	!!! Processing in batch mode. !!!
C	Version 1.2.3, April 28, 1997 by M. Knappertsbusch,
C	modified to batch mode on 21.12.1999
C	modified to write errorneous images to external file. 21.1.2000
C	modified to record missing images, that are read from FILE_LIST. 24.1.2000
C	modified by changing MAG from INTEGER to REAL. 18.12.2000
C	modified by changing the calibration (New XPREC and YPREC). 22.12.2000
C	modified from Trace33_batch.f for application with AMOR using Sony DXC-390P camera. 30.5.2007
C	Modified from Trace_AMOR1_batch.f to Trace_AMOR2_batch.f: Correction of wrong YPREC. 13.03.2018
C
C	Program to trace the outline of a series of digitized light images
C	in batch mode. Input images are grey-level pictures with 480 lines
C	and 680 pixels per line.
C
C
C	!!!!! Output are the X,Y coordinates of the outline
C	in COUNTER-CLOCKWISE direction !!!!!
C	Program includes conversion of pixels in micrometers.
C	Program calibrated to the Leica MZ6 binocular stereo microscope.
C
C	Input files are:
C		- A File with name FILE_LIST containing a list of all
C		  names of the input images (filename) to be traced,
C		  and the magnifications, that were applied to each image.
C		  FILE_LIST is a character variable and can bear any name
C		  up to 20 characters in length.
C
C		- All image files to be traced. The names of these raw images
C		  are interpreted through the character variable INPUT.
C		  INPUT is exactly 5 characters long.
C
C		- The program calls for the sample name through the
C		  character variable SAMPLE. SAMPLE is esactly 11 characters
C		  long (including orientation codes K or U for Keel or Umbilical).
C
C		  INPUT and SAMPLE are used to compose the name of the corresponding
C		  output file (suffix _T).
C
C	Output files are:
C		- One file OUTPUT for every traced sample, containing the x,y 
C		  coordinates of the outline. The character variable OUTPUT is
C		  17 characters long (including the suffix _T).
C
C	Note: All file handling for morphometric analyses follows a name convention:
C
C	filename      = Name of sample, 15 characters long.
C	filename_T    = Input file for raw curves with cartesian coordinates, 17 characters long.
C	filename_POL  = Polar coordinates of filename,                        19 characters long.
C	filename_INT  = Interpolated (cartesian) coordinates of filename,     19 characters long.
C	filename_NORM = Normalized curve,                                     20 characters long.
C
C
C***	Note: MAG must be real, otherwise results may be wrong !
C
	INTEGER I,N,N0,X,Y,THRES,NPTS,NF
	INTEGER GREY(1:8), K(1:8), A(1:8,1:8)
	LOGICAL IMAGEX
	REAL MAG
	DOUBLE PRECISION XPREC,YPREC,XRES,YRES
	CHARACTER*1 CHAR(1:8),ANSWER,COMMA
	CHARACTER*5 INPUT
	CHARACTER*17 OUTPUT
	CHARACTER*11 SAMPLE
	CHARACTER*20 FILE_LIST
	PARAMETER (THRES=49)
C
C	In A(I,J) is the local numbering of points:
C
	DATA ((A(I,J),I=1,8),J=1,8)/1,
     *    2,3,4,5,6,7,8,
     *  2,3,4,5,6,7,8,1,
     *  3,4,5,6,7,8,1,2,
     *  4,5,6,7,8,1,2,3,
     *  5,6,7,8,1,2,3,4,
     *  6,7,8,1,2,3,4,5,
     *  7,8,1,2,3,4,5,6,
     *  8,1,2,3,4,5,6,7/
C 
C	Experimental setup and calibration:
C
C       A Leica MZ6 Stereo-Binocular eqipped with a Sony DXC-390P colour
C	camera and using the AMOR imaging system (PC based). The camera is connected to the
C	microscope via a 1x Cmount, and Achromat 1x and 2x objectives.
C	The magnification to be entered is calculated as the product
C	of the objective lens inserted (1x or 2x) times the number indicated
C	at the zoom-knob. Example: if the objective lens is 1x and the
C	number at the zoom-knob indicates 4.0 then the magnification to be
C	entered is 4.0.
C
C	In the following factors are calculated to convert pixels into µm with the 
C	CCD camera, using the NTSC (640x480 pixels format in NIH):
C	(XPREC and XRES are in horizontal, YPREC and YRES are in vertical direction).
C	This calibration includes now a newer one, which was done with the 0.01mm scale
C	(=small scale) and the 0.1mm scale (large scale from Leica). This has not much
C	changes in the final result but the linearity of XPREC and YPREC over the entire
C	range of magnifications from 0.63x to 8x could so be demonstrated (see Knappertsbusch
C	2004-2007, MorphCol supplement #4, @ conclusion 7.
C	Equation (II) was wrongly implemeted in Trace_AMOR1_batch.out.
C	See especially note "Fehlerhaftes Trace_AMOR1_batch.out Programm", in docu
C	"Weitergehende PlŠne fŸr die Entwicklung von AMOR" from 13 March 2018.
C
C	XPREC = 0.12986 * MAG + 0.00022300    (in pixel/µm), equation (I)
C	YPREC = 0.15742 * MAG - 0.00023756    (in pixel/µm), equation (II)
C
C	The necessary resolutions XRES and YRES are the calculated as
C
C	XRES = 1/XPREC     (in µm/pixel)
C	YRES = 1/YPREC     (in µm/pixel)
C
C	Specifications on computer side (PC with AMOR program):
C	Frame grabber: IMAQ 1405 from National Instruments.
C	Video control: Grey-scale, full-size, NTSC format
C	Digital image: resolution: 72 pixel per inch, image size: 640x480 pixel
C	
C
C
	WRITE(9,*) '               Program Trace'
	WRITE(9,*) 'Trace_AMOR2_batch.out, calibration from 11.04.2007'
	WRITE(9,*) '            Processing in batch mode'
	WRITE(9,*) ' '
	WRITE(9,*) '. . .Enter sample name (11 chars, with K or U). . .'
	READ(9,102) SAMPLE
 102	FORMAT(A11)
C
	WRITE(9,*) '. . .Output in Pixels (1) or Microns (2) ?. . .'
	READ(9,4) ANSWER
 4	FORMAT(A1)
 	IF (ANSWER.EQ.'1') THEN
	  WRITE(9,*) 'Output is in pixels'
	ELSE
	  WRITE(9,*) 'Output is in µm'
	  WRITE(9,*) ' '
	END IF
C
C	Get filenames of images and magnifications from FILE_LIST, and start opening
C	and reading for first image:
C
	WRITE(9,*) '. . .Enter list containing image
     1 names (20 chars) and MAG. . .'
     	WRITE(9,*) 'Format: 5 chars, comma ,1 real value (2 decimals)'
	WRITE(9,*) 'Note: MAG=Objective lens x zoom factor'
	WRITE(9,*) 'Example: 0101r,3.20 (specimen=0101r,
     2 magnification=3.20x)'
	READ(9,101) FILE_LIST
 101	FORMAT(A20)
C
C
C
C
C	Initialization of NF (NF= number of files traced):
C
	NF=0
C
C	Recording errors during program run:
C
	OPEN(21,FILE='Image_Errors',STATUS='NEW')
	WRITE(21,96) SAMPLE
 96	FORMAT('Errors in Images of Sample: ',A11)
 	OPEN(17,FILE='Missing_images',STATUS='NEW')
C
C	Now reading filenames etc:
C
	OPEN(20,FILE=FILE_LIST,STATUS='OLD')
 100	READ(20,110,END=999) INPUT,COMMA,MAG
 110	FORMAT(A5,A1,F4.2)
C
C	Check first, whether image with name INPUT exists or not. If it does not exist:
C	print a message and continue reading next image.
C
	INQUIRE(FILE=INPUT,EXIST=IMAGEX)
	IF (.NOT.IMAGEX) THEN
 	  WRITE(9,*) 'Error: image ',INPUT,' is missing'
	  WRITE(17,97) INPUT
 97	  FORMAT(A5)
	  WRITE(9,*) '. . .Program continues. . .'
 	  GOTO 100
	END IF
C
C	Counting the number of images (=NF) treated:
C
	NF=NF+1
C
C	Initialization of NPTS: (=Number of points in outline):
C
	NPTS=1
C
C	Determination of the names of the OUTPUT files (17 chars long), and
C	writing information to screen:
C
	OUTPUT=(SAMPLE(1:11)//INPUT(1:4))//'_T'
	WRITE(9,111) OUTPUT,MAG
 111	FORMAT('  Specimen: ',A17,' at magnification: ',F4.2,'x')
C
C
C	Determination of XRES and YRES (has been calibrated using the NTSC 
C	(640x480 pixel) format in NIH-Image:
C
	XPREC = 0.12986*MAG+0.000223
	YPREC = 0.15742*MAG-0.00023756
C
	XRES=1/XPREC
	YRES=1/YPREC
C
C
C	WRITE(9,*) '. . .calculating. . .'
C
	OPEN(15,FILE=INPUT,ACCESS='DIRECT',
     *  RECL=1,FORM='FORMATTED',STATUS='OLD')
     	OPEN(16,FILE=OUTPUT,STATUS='NEW')
C
C	Find first point of outline in the middle of the
C	TV-screen (this is line 240; its record is calculated
C	as 152960=(240-1)*640):
C
	DO 1, I=1,640
	  N=152960+I
	  READ(15,50,REC=N) CHAR(1)
	  GREY(1)=ICHAR(CHAR(1))
	    IF (GREY(1).GE.THRES) THEN
	      N0=N
	      GOTO 5
	    END IF
 1	CONTINUE
C
 5	CALL COORD(N,X,Y)
C 	WRITE(9,*) X,Y
C
C	Output of results to file:
C
	IF (ANSWER.EQ.'1') THEN
C	  Output of X and Y in pixels:
 	  WRITE(16,*) X,',',Y
	ELSE
C	  Output of X and Y in micrometers:
	  WRITE(16,*) X*XRES,',',Y*YRES
	END IF
C
C
C	STOP Criterium if there is an error in the image:
C	(NPTS indicates the number of points written to the outline file; usually,
C	an outline of forams at a magnification of ²4x has less than 2500 points).
C	The errorneous INPUT images are written into external file with name
C	"Image_Errors", that can be used as input for subsequent runs.
C
	NPTS=NPTS+1
	IF (NPTS.GT.2500) THEN
	  WRITE(9,119) OUTPUT	 
	  WRITE(21,110) INPUT,COMMA,MAG
 119	  FORMAT('Error in image ',A17,'. Reading next image')
	  GOTO 100
	END IF
C
C
C
C
C
C
C	Local initialization: At each point (pixel) of the outline (record number N)
C	the eight surrounding pixels are read, then their grey-levels
C	are determined, and the pixels are numbered from 1 to 8:
C
C
C		2	3	4
C
C		1	N	5
C
C		8	7	6
C
C
	READ(15,50,REC=N-1)   CHAR(1)
	READ(15,50,REC=N-641) CHAR(2)
	READ(15,50,REC=N-640) CHAR(3)
	READ(15,50,REC=N-639) CHAR(4)
	READ(15,50,REC=N+1)   CHAR(5)
	READ(15,50,REC=N+641) CHAR(6)
	READ(15,50,REC=N+640) CHAR(7)
	READ(15,50,REC=N+639) CHAR(8)
 50	FORMAT(A1)
C
C	Determination of the grey-levels:
C
	DO 10, I=1,8
	  GREY(I)=ICHAR(CHAR(I))
 10	CONTINUE
C
C	Assignment of the local pixel number (i.e. 1 through 8)
C	to the record number N for each pixel.
C	(i.e. in K are the record numers of the pixels 1 to 8 stored,
C	with K having the idices from i to 8):
C
	K(1)=N-1
	K(2)=N-641
	K(3)=N-640
	K(4)=N-639
	K(5)=N+1
	K(6)=N+641
	K(7)=N+640
	K(8)=N+639
C
C	Search algoritm for the next point of the outline: The algoritm
C	searches first for the local point to start with: This first
C	point has the grey-level ³ Threshold. From this initial
C	local point it searches the next point, where the greylevel is
C	lower than Threshold. The pixel, which is previous to that
C	point is the next point of the outline (i.e. the next center for
C	the search algoritm).
C
C	Look for first point of lokal environment-search:
C
	DO 20, I=1,8
	  IF (GREY(I).GE.THRES) THEN
C
C	    First point of search path found. Check every next surrounding point in
C	    clockwise direction:
C
	    DO 30, J=1,8
	      IF (GREY(A(I,J)).LT.THRES) THEN
	        IF (A(I,J).EQ.1) THEN
		  N=K(8)
		ELSE
		  N=K(A(I,J)-1)
		END IF
C
C	          Stop-Criterium:
C
		  IF (N.EQ.N0) THEN
C		    WRITE(9,*) 'Outline closed'
		    CLOSE(16)
		    CLOSE(15)
C
C		    Reading next file:
		    GOTO 100
C
		  END IF
C
		  GOTO 5
		END IF
 30	    CONTINUE
	  END IF
 20	CONTINUE
C
 999	WRITE(9,*) ' '
 	WRITE(9,*) NF,' images traced and written to files'
	WRITE(9,*) 'One file "Image_Errors" written'
	WRITE(9,*) 'One file "Missing_images" written'
 	WRITE(9,*) 'Enter a key to end the program'
 	PAUSE 888
 888	CONTINUE
C
 	STOP
	END




	SUBROUTINE COORD(M,U,V)
C
C	Determines the cartesian coordinates U and V as a function of
C	the record number M representing a pixel in the image.
C	The image has 480 lines and 640 pixels per line.
C
C	The record number starts with M=1 (representing the first
C	pixel in line 1).
C
	INTEGER M,U,V
C
C	This determination of U and V is valid for all pixels except those
C	that are at the end of each line (i.e. pixel no. 640, 1280,...).
C
	U=MOD(M,640)
	V=479-INT(M/640)
C
	RETURN
	END