	program peakmatch
C first, read everyone into memory
C  April 17, 2002
C Phil Massey, Lowell Observatory
C compares the size of the peak stellar intensities and determines relative
C scale factors for use by imcombine
	implicit double precision (a-h,o-z)
	parameter(maxstars=250000)
	parameter(maxfiles=15)
	dimension ra(maxstars,maxfiles),dec(maxstars,maxfiles),
     &            peak(maxstars,maxfiles),nstars(maxfiles)
	character*80 name,inline
	character*60 names(maxfiles)
	dimension fact(maxstars)

*  SECONDS TO RADIANS
      PARAMETER (S2R=0.72722052166430398D-4)
*  ARC SECONDS TO RADIANS
      PARAMETER (AS2R=0.48481368110953599D-5)
C
*  SECONDS TO RADIANS
      STOR(RH,RM,RS)=S2R*(60D0*(60D0*RH+RM)+RS)
*  ARC SECONDS TO RADIANS
      ASTOR(DD,DM,DS)=AS2R*(60D0*(60D0*DD+DM)+DS)
	data hit/0.7/

C Read everyone into memory
C First, open the list of image names:
	open(unit=1,file='temppeaks',status='old')
	ik=0
1	read(1,2,end=99) name
2	format(a80)
	last=lnblnk(name)
	open(unit=2,file=name(1:last)//'prd',status='old')
	ik=ik+1
	if(ik.gt.maxfiles) stop "too many files"
	names(ik)=''
	names(ik)=name(1:last)
	istar=0
3	read(2,2,end=88) inline
	call bear(inline,ihr,imin,sec,ideg,iamin,asec,isuc)
	if(isuc.eq.0) go to 3
	istar=istar+1
	if(istar.gt.maxstars) stop "too many stars"
	ra(istar,ik)=stor(ihr,imin,sec)
	dec(istar,ik)=astor(ideg,iamin,asec)
	ic=0
	do 4 i=1,80
	if(inline(i:i).eq.':') ic=ic+1
	if(ic.lt.3) go to 4
	go to 5
4	continue
	stop 'problem here'
5	read(inline(i+9:80),*) peak(istar,ik)
	go to 3
88	close(unit=2)
	nstars(ik)=istar
	write(6,89) ik,nstars(ik),names(ik)
89	format('For file',i3,' found', i7, ' stars ',a60)
	go to 1
99	nfiles=ik
C Now, for every star in 1, try to find a match in the other files
C Now for each file:
	factor=1.000
	write(3,207) factor,names(1)
	do 200 j=2,nfiles
	ifact=0
	do 205 l=1,nstars(j)
	do 100 i=1,nstars(1)
	call apa(ra(i,1),dec(i,1),ra(l,j),dec(l,j),sep,pa)
	dist=sep/as2r
	if(dist.lt.hit) then
		ifact=ifact+1
		fact(ifact)=peak(i,1)/peak(l,j)
		write(14,*) fact(ifact)
		go to 205
	end if
100	continue
205	continue
	call median(fact,ifact,amed)
	write(6,206) j,ifact,amed
206	format('Median factor between file 1 and file',2i5,
     &  f12.4)
	write(3,207) amed,names(j)
207	format(f7.4,1x,a60)
200	continue
	stop
	end
        subroutine apa(rac,dcc,rao,dco,dist,pa)
c computes the distance (delt) and position (pa) along the
c sky FROM (rac,dec) to (rao,dco)
c written aug 14, 1987
C rewritten August 21, 1989, to do using "standard coords", rather
C than the old and peculiar (but possibly equivalent!) way
c plm
      implicit double precision (a-h,o-z)
        parameter (PI=3.141592653589793238462643)
C CONVERT TO PROPER COORDS
        CALL DS2TP(RAO,DCO,RAC,DCC,XI,ETA)
        pa=0.
        if(xi*eta.ne.0.)PA=DATAN2(XI,ETA)
        if(pa.le.0.) pa=pa+2.*pi
        DIST=DSQRT(XI**2+ETA**2)
        return
        end
        include '/tofu/home/massey/observing/astromsub.f'
        include 'mediandbl.f'
	include 'bearsub.f'
	
