C+
C NAME:
C	A_mainf_transit_test
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 A_mainf_transit_test
	implicit none

	! Dimensions for node, nodeid, hits above must be at least 2*(4**(level+1))
	! For LEVEL = 12, KMAX = 134217728
	! For LEVEL =  8, KMAX =    524288

	integer*2	  TOP_LEVEL
	parameter	( TOP_LEVEL = 10 )

	integer*8	  TOP_KMAX
	parameter	( TOP_KMAX = 2*(4**(TOP_LEVEL+1)) )

	integer*2	level			/TOP_LEVEL/
	integer*8	kmax			/TOP_KMAX/

	real*4		frame(1272,256)
C	real*4		frame(636,128)
C	real*4		frame(318,64)
	real*4  	ped			/0./
	real*4		dark			/0./
	real*4		hroom			/49400./		!'headroom'=new onboard value of 65k ADUs
	real*4		adu2e
	real*8		affine(6)		/6*0.0D0/		! Affine parameters zero for now
	real*8		cra
	real*8		cdec
	real*8		alfa
	real*8		beta
	real*8		gamma
	real*8		euler(3)
	real*8		camlat
	real*8		degperminute
	real*8		ra0
	real*8		arg
	real*8		arg1
	real*8		dusk

	integer*4	mode			/1/			!onboard binning mode
	integer*4	ic			/3/			!Camera #
	integer*4	ifm			/0/			!ROI mask (0=no,1=yes)
	integer*4	id			/1272/
	integer*4	jd			/256/
	integer*4	icrcnt			/0/
	integer*4	neighbors		/0/
	integer*4	ihist(250)
	integer*4	i
	integer*4	j

	character	outfile*37,listfile*23,framename*25,infile*45,starname*17
c
c	AARON'S DECLARATIONS
c
	integer*4 idarkcr, dumcnt4, icra, goodframe, ss, date, maxk
	integer*4 minf, maxf, switcher
	integer*4 k, kk, loopcnt, starcount, dumcount
	real*8 response, arg, arg1, phi, rlat
	real*8 qcam(3,4), qc(4), qsc(4), qandy(4), qdum(4), qfinal(4)

c	date=0225
	date=0301
	write(listfile,'(A15,I4.4,A4)') 'DATA_SMEI/list_',date,'.txt'
	!!Camera 1
	if(date.eq.0225)then
	  maxk=6629
	  write(outfile,'(A37)') 'Star_Transits/0225_01_starlist_01.txt'
	endif
	!!Camera 3
	if(date.eq.0301)then
	  maxk=3271
	  write(outfile,'(A37)') 'Star_Transits/0301_01_starlist_01.txt'
	endif
	open(11,file=outfile,status='new')

	!! loop through all of the stars
	open(68,file='DATA_SMEI/BrightStarsSorteda.txt',status='old')
	starcount=0
	do ss=1,53
	read(68,69)cra,cdec,starname
69	format(f12.8,2x,f12.8,49x,A17)
	if(cdec.gt.60.0.or.cdec.lt.-60.0) goto 976
	
	!!Get rid of crap in the name
	do k=1,17
	  if(starname(k:k).eq.' '.or.starname(k:k).eq.'*') write(starname(k:k),'(A1)') '_'
	enddo

	print*, starname, cra, cdec
	dumcount=0

	!! Quaternion from Andy's Euler angles for camera 3
	qcam(3,1) = 0.216504D0
	qcam(3,2) = 0.432944D0
	qcam(3,3) = -0.868351D0
	qcam(3,4) = -0.107945D0
	!! refined camera 1 quaternions, ap003, a3tx9, tx995
	qcam(1,1) =-0.2377267D0
	qcam(1,2) = 0.4680082D0
	qcam(1,3) = 0.8481501D0
	qcam(1,4) = 0.0713853D0
	!! updated test quaternions(gp020,bm0035) with x0 shift, thetax & thetay correction
	qcam(2,1) =-0.26839441D0
	qcam(2,2) = 0.01748647D0
	qcam(2,3) = 0.96205384D0
	qcam(2,4) = 0.04594647D0

	!!quaternion to rotate Andy's coordinates to Don's
	qandy(1)=0.5D0
	qandy(2)=0.5D0
	qandy(3)=0.5D0
	qandy(4)=0.5D0
	data ihist /250*0/

	open(13,file=listfile,status='old')

	loopcnt=0
	switcher=0
	do k=1,maxk

	  read(13,13)framename
13	  format(a25)

	  if(k.gt.0)then
	  goodframe=0

	        write(infile,'(A15,I4.4,A1,A25)') 'DATA_SMEI/2003_',date,'/',framename
                if(mod(k,500).eq.0) print*,k

C		call read_frame(infile,mode,ic,id,jd,frame,qdum)
		call read_framecrx_old_i2(infile,mode,ic,id,jd,frame,qdum)

		!! input quaternion is of form ai+bj+ck+d, not d+ai+bj+cK -> swap with new quaternion
		qsc(1)=qdum(4)
		qsc(2)=qdum(1)
		qsc(3)=qdum(2)
		qsc(4)=qdum(3)

		!! multiply camera quaternion into spacecraft(frame) quaternion and take the inverse
		!! of camera quaternion because it goes from sky to ccd
		qc(1)=qcam(ic,1)
		qc(2)=-qcam(ic,2)
		qc(3)=-qcam(ic,3)
		qc(4)=-qcam(ic,4)

		call qmultiply(qsc,qc,qdum)
		call qmultiply(qdum,qandy,qfinal)
		
		!!take inverse of final quaternion: going from ccd to sky, not the other way around
C		qfinal(1)=qfinal(1)
		qfinal(2)=-qfinal(2)
		qfinal(3)=-qfinal(3)
		qfinal(4)=-qfinal(4)

C		  find ra/dec of center of FOV for a frame
C		  find ra of center of the first frame; use to calculate ra bounds in final image
C		  car/cdec calculated from thetax=thetay=0
C	              phi=0.D0
C		      rlat=90.D0
C	              call rotateq(qfinal,phi,rlat) !!rotate to standard system
C		      print*, '		RA/DEC:  ',phi,rlat

		call indexframe(mode,ic,id,jd,affine,qfinal,goodframe,cra,cdec)
		if(goodframe.eq.1) then
		  print*, switcher,k,': ',infile
		  switcher=switcher+1
		  if(switcher.eq.1)  minf=k
		  maxf=minf+switcher-1
		else
		  if(switcher.gt.5)then
		    write(11,23) cra,cdec,minf-3,maxf+3,starname
23		    format(f10.6,2x,f10.6,2x,i5,2x,i5,2x,A17)
		    dumcount=dumcount+1
		  endif
		  switcher=0
		endif

	  endif
	  enddo
	  if(dumcount.gt.0) starcount=starcount+1
	  if(starcount.eq.15) goto 99

976	  close(13)
	enddo
99	close(11)

	print*, starcount, 'Stars Tracked!'
	print*, 'Fortran out!'

	end