C+
C NAME:
C	transit_centroider
C PURPOSE:
C	to find the centroid of a star averaged over a transit
C CATEGORY:
C	SMEI Data processing
C INPUTS:
C
C OUTPUTS:
C
C CALLS:
C       none
C PROCEDURE:
C
C MODIFICATION HISTORY:
C	April, 2003 A. Smith (UCSD/CASS)
C-
	program transit_centroider
	implicit none

	character	infile*17, outfile*16, dum*20
	integer*4	image(100,100), i, j, x0, y0, maxr,minr,response, test
	integer*8	xsum, ysum, sum
	real*8		xcen, ycen, ra, dec,dra, ddec, cra, cdec
	real*4		radius
	
	print*, 'Fortran in...'

	xsum=0
	ysum=0
	sum=0
	maxr=0
	minr=60000

CCC	Coordinates and Files for Aldebaran
	cra=68.9802D0
	cdec=16.5093D0
	write (infile,'(A17)') 'Aldebar_tran1.grd'
	write (outfile,'(A16)') 'Aldebar_cent.grd'
CCC	Coordinates and Files for HD 20902
c	cra=51.0807D0
c	cdec=49.8611D0
c	write (infile,'(A17)') 'HD20902_tran1.grd'
c	write (outfile,'(A16)') 'HD20902_cent.grd'
CCC	Coordinates and Files for Sirius (HD 48915)
c	cra=101.28717D0
c	cdec=-16.71611D0
c	write (infile,'(A17)') 'Sirius__tran1.grd'
c	write (outfile,'(A16)') 'Sirius__cent.grd'
CCC	coordinates and files for Vega (alfa lyra)
c	cra=279.23475D0
c	cdec=38.78361D0
c	write (infile,'(A17)') 'Vega____tran1.grd'
c	write (outfile,'(A16)') 'Vega____cent.grd'
CCC	coordinates and files for Alcaid
c	cra=206.8852D0
c	cdec=49.3133D0
c	write (infile,'(A17)') 'Alcaid__tran1.grd'
	write (outfile,'(A16)') 'Alcaid__cent.grd'
CCC	coordinates and files for Eta Draconis
c	cra=245.99786D0
c	cdec=61.51421D0
c	write (infile,'(A17)') 'EtaDrac_tran1.grd'
c	write (outfile,'(A16)') 'EtaDrac_cent.grd'
CCC	coordinates of Polaris
c	cra=37.752D0
c	cdec=89.264D0
c	write (infile,'(A17)') 'Polaris_tran1.grd'
c	write (outfile,'(A16)') 'Polaris_cent.grd'


	!! read in headers, then image array
	open(12,file=infile,form='formatted',readonly)
	read(12,120)dum
	read(12,120)dum
	read(12,120)dum
	read(12,120)dum
	read(12,120)dum
120	format(a20)
	read(12,12)((image(i,j),i=1,100),j=1,100)
12	format(20i6)
	close(12)

	!! find brightest pixel to draw box around
	do j=1,100
	 do i=1,100
	   if(image(i,j).gt.maxr)then
	     maxr=image(i,j)
	     x0=i
	     y0=j
	   endif
	 enddo
	enddo
	
	y0=y0-3
	
	!!find background within circle
	do j=1,100
	 do i=1,100
	   radius = sqrt( float((i-x0)*(i-x0) + (j-y0)*(j-y0)) )
	   if(radius.lt.18)then
	     if(image(i,j).lt.minr) minr=image(i,j)
	   endif
	 enddo
	enddo
	
	!! find centroid of star
	do j=1,100
	 do i=1,100
	   radius = sqrt( float((i-x0)*(i-x0) + (j-y0)*(j-y0)) )
	   if(radius.lt.18)then
	      response=image(i,j)-minr
	      sum = sum + response
	      xsum = xsum + (i*response)
	      ysum = ysum + (j*response)
	   endif
	 enddo
	enddo
	xcen=dfloat(xsum)/dfloat(sum)
	ycen=dfloat(ysum)/dfloat(sum)

	image(nint(xcen),nint(ycen))=-100

	ra = cra - (xcen - 50.5D0)*0.5D0
	dec = cdec + (ycen - 50.5D0)*0.5D0

	ddec=dec-cdec
	dra=ra-cra
	
	print*, 'X,Y:     ',xcen,ycen
	print*, 'RA/DEC:  ',ra,dec
	print*, 'Discrep: ',dra,ddec
	
	!!Write punctured image to file
	open (10,file = outfile, form = 'formatted')
	write(10,10)
10	format('DSAA'/'100 100'/'0 100'/'0 100'/'0 10000')
	do j=1,100
	  write(10,'(20i6)')(image(i,j),i=1,100)
	enddo
	close(10)

	print*, 'Fortran out!!!'

	end