      program ukscat_www
c+
c f77 ukscat_www.f -o ukscat_www.e
c           this is a rather crude (but working) program
c           which will search the uk schmidt plate catalogue
c           and select plates according to the following
c           criteria.
c                 ra and dec of plate centre
c                 plate number
c                 waveband
c                 exposure time
c                 prism or direct plates
c                 non-survey reference code ('t number')
c                 field number
c                 appearance of specified ra and dec on plate
c
c           it must be run via the command file ukscat.com.
c           it is interactive and prompting, and is designed
c           to be reasonably easy to run from a remote node.
c         
c           by k p tritton and m j vickerstaff
c           last modified 22 feb 1983
c           modified 1996 aug 9 -mjd for command-line arguments and WWW
c-


      character*1 order,band,prism,sband,sprism
      character line*70,str1*1,str3*3,numfs*5,str22*22,argv*40
      character*4  s7,s10
      character*5  s3,s4,s5
      character*6  s6,s8,s9
      character*7  s2
      character*8  s1,s11

      integer lora,hira,lodec,hidec,out1,out2,ptra,ptdec
      double precision radiff,sra,sdec,sdeca,pra,pdec,pdeca
      double precision pi,xi,eta,xpl,ypl
      logical debug

      debug=.false.
c      debug=.true.
      c=3073.0
      pi=3.1415926
      in1=5
      in2=7
      out1=6
      out2=6
c      open(unit=out2,file='search.lis')

c           explain options
c
c     write (6,*) 'enter search option:'
c     write (6,*) '(0) specify order (Plate no/Ra) [P]
c     write (6,*) '(a) plate number  (min max)            [0,0]
c     write (6,*) '(b) plate center (rmin,rmax,dmin,dmax) [0,0,0,0]
c     write (6,*) '(c) waveband (BUV... X=other 0=any)    [0]
c     write (6,*) '(d) exposure time (min,max)            [0,0]
c     write (6,*) '(e) prism or direct (Prism/Dir/0=any)  [0]
c     write (6,*) '(f) t number        (0=any)            [0]
c     write (6,*) '(g) field number    (0=any)            [0]
c     write (6,*) '(h) plates containing given ra, dec    [-1,0]

      order = 'P'
      npmin = 0
      npmax = 0
      iramin = 0
      iramax = 0
      band = '0'
      nexpmin = 0
      nexpmax = 0
      prism = '0'
      ntnum = 0
      nfield = 0
      ptra = -1

c *** try for presence of command line arguments
      m=iargc()
      ifu=0

c *** if reasonable number of arguments try and interpret them
      if(m.eq.9)then
        call getarg(1,argv)
        read(argv,'(a)',err=4,end=4) order
        call getarg(2,argv)
        read(argv,*,iostat=ioerr) npmin,npmax
        if (ioerr.ne.0) then
          npmin=0
          npmax=0
        end if
        call getarg(3,argv)
        read(argv,*,iostat=ioerr) iramin,iramax,idemin,idemax
        if (ioerr.ne.0) then
          iramin=0
          iramax=0
          idemin=0
          idemax=0
        end if
        call getarg(4,argv)
        read(argv,'(a)',err=4,end=4) band
        call getarg(5,argv) 
        read(argv,*,iostat=ioerr) nexpmin,nexpmax
        if (ioerr.ne.0) then
          nexpmin=0
          nexpmax=0
        end if
        call getarg(6,argv)
        read(argv,'(a)',err=4,end=4) prism
        call getarg(7,argv)
        read(argv,*,iostat=ioerr) ntnum
        if (ioerr.ne.0) then
          ntnum=0
        end if
        call getarg(8,argv)
        read(argv,*,iostat=ioerr) nfield
        if (ioerr.ne.0) then
          nfield=0
        end if
        call getarg(9,argv)
        read(argv,*,iostat=ioerr) ptra,ptdec
        if (ioerr.ne.0) then
          ptra=-1
          ptdec=0
        end if
        goto 5
      endif
 4    ifu=1
 5    continue
      if(m.ne.9.or.ifu.eq.1) call all_crap( 'parameters wrong?' )
      if (idemax.lt.idemin) call all_crap( 'dec maximum less than 
     +minimum')
      if (npmax.lt.npmin) call all_crap( 'pnum maximum less 
     +than minimum')
      if (nexpmax.lt.nexpmin) call all_crap( 'exposure max less 
     +than min')
      if ((ptra.gt.24000).or.(ptdec.lt.-9000).or.(ptdec.gt.9000))
     +   call all_crap('invalid point co-ordinates')

      if (debug) then
         write (6,*) 'Parameters    ',m
         write (6,*) 'Order         ',order
         write (6,*) 'Plate numbers ',npmin,npmax
         write (6,*) 'Plate centres ',iramin,iramax,idemin,idemax
         write (6,*) 'Band          ',band
         write (6,*) 'Exposure      ',nexpmin,nexpmax
         write (6,*) 'Prism/direct  ',prism
         write (6,*) 'T number      ',ntnum
         write (6,*) 'Field number  ',nfield
         write (6,*) 'Position      ',ptra,ptdec
      end if

c *** open the catalogue file    *************************
      ioerr=0
c      open(in2,file='ukst_catalogue',status='old',iostat=ioerr)
      open(in2,file='/home/mar/sss/catlog.lis'
     +  ,status='old',iostat=ioerr)
      if (ioerr.ne.0) call all_crap( 'cant open catalogue file')
      read(in2,'(a)',iostat=ioerr) line
      do while (ioerr.eq.0)
        read(in2,'(a)',iostat=ioerr) line
      end do
      write(out2,'(x,''<HTML> <TITLE>'')')
      write(out2,'(x,''UKST Plate Catalogue Search Results'')')
      write(out2,'(x,''</TITLE>'')')
      write(out2,'(x,''<H1><B>'')')
      write(out2,'(x,''UKST Plate Catalogue Search Results'')')
      write(out2,'(x,''</B></H1>'')')
      write(out2,'(x,''<P>'')')
      write(out2,'(x,''Please go'')')
      write(out2,'(x,''<a 
     >href="http://www.roe.ac.uk/ukstu/pltcat.html">'')')
      write(out2,'(x,''here</a> for a description '')')
      write(out2,'(x,''of the catalogue.'')')
      write(out2,'(x,''<H3><B>'')')
      write(out2,'(x,''Complete to plate '',a'' taken ''a''/''a''/''a)') 
     +     line(1:8),line(35:36),line(33:34),line(31:32)
      write(out2,'(x,''</B></H1>'')')
      write(out2,'(x,''<P>'')')
      write(out2,'(x,''<PRE>'')')
      write(out2,'(x,''Order         ''a)')order
      write(out2,'(x,''Plate numbers ''2i6)')npmin,npmax
      write(out2,'(x,''Plate centres ''i6.5,i6.4,i6.5,i6.4)')
     +     iramin,iramax,idemin,idemax
      write(out2,'(x,''Band          ''a)')band
      write(out2,'(x,''Exposure      ''2i6)')nexpmin,nexpmax
      write(out2,'(x,''Prism/direct  ''a)')prism
      write(out2,'(x,''T number      ''i6)')ntnum
      write(out2,'(x,''Field number  ''i6)')nfield
      write(out2,'(x,''Position      ''i6.5,i6.4)')ptra,ptdec
      write(out2,'(x)')
ccccc      write(out2,'(x,'' Plate    Programme or    RA   dec    Date   LST''
ccccc     +     '' Emulsion   Exposure  Grade'')')
      write(out2,'(x,'' Plate    T no. or Field  RA   
     +dec    Date   LST''
     +     ''  Emulsion   Exposure  Grade'')')
      write(out2,'(x,'' number   survey              
     +                  ''
     +     ''  and filter     time'')')
      write(out2,'(x)')

      close(in2)

      ioerr=0
      if (order.eq.'R'.or.order.eq.'r') then
c      open(in2,file='ukst_ra_catalogue',status='old',iostat=ioerr)
      open(in2,file='/home/mar/sss/catra.lis',iostat=ioerr)
      else
c      open(in2,file='ukst_catalogue',status='old',iostat=ioerr)
      open(in2,file='/home/mar/sss/catlog.lis',iostat=ioerr)
      end if
      if (ioerr.ne.0) call all_crap( 'cant open catalogue file')

c *** start search
      icount=0
      nread=0

c *** read the catalogue file in *************************
 204  read (in2,'(a)',end=999,err=998) line
      ifail=0
      nread=nread+1

      read (line,'(a1,a1,i5,a1,a3,i4,a5,2(i5),a22,i4)',err=998)
     + str1,sband,num,sprism,str3,numt,numfs,ira,idec,str22,numx

      if (npmin.ne.0.or.npmax.ne.0) then
        if (num.lt.npmin.or.num.gt.npmax) ifail = 1
      end if

      if (iramax.ne.0.or.iramin.ne.0) then
      if (iramax.ge.iramin) then
        if (ira.lt.iramin.or.ira.gt.iramax.or.idec.lt.idemin.or.idec.
     +    gt.idemax) ifail = 1
      else
        if ((ira.lt.iramin.and.ira.gt.iramax).or.idec.lt.idemin.or.
     +    idec.gt.idemax) ifail = 1
      end if
      end if

      if (band.ne.'0') then
      if (band.eq.'u') band = 'U'
      if (band.eq.'b') band = 'B'
      if (band.eq.'v') band = 'V'
      if (band.eq.'r') band = 'R'
      if (band.eq.'z') band = 'Z'
      if (band.eq.'j') band = 'J'
      if (band.eq.'i') band = 'I'
      if (band.eq.'x') band = 'X'
      if ((band.eq.'U'.or.band.eq.'B'.or.band.eq.'V'.or.band.eq.'R'.
     +or.band.eq.'Z'.or.band.eq.'J'.or.band.eq.'I').and.sband.ne.band)
     +ifail=1
      if (band.eq.'X'.and.((sband.eq.'U'.or.sband.eq.'B'.or.
     +  sband.eq.'V'.or.sband.eq.'R'
     +  .or.sband.eq.'Z'.or.sband.eq.'J'.or.sband.eq.'I'))) ifail = 1
      end if

      if (nexpmin.ne.0.or.nexpmax.ne.0) then
        if (numx.lt.nexpmin.or.numx.gt.nexpmax) ifail = 1
      end if

      if ((prism.eq.'P'.or.prism.eq.'p').and.sprism.ne.'P') ifail = 1
      if ((prism.eq.'D'.or.prism.eq.'d').and.sprism.eq.'P') ifail = 1

      if (ntnum.ne.0.and.ntnum.ne.numt) ifail = 1

      if (nfield.ne.0) then
        read(numfs,*,iostat=ioerr) numf
        if (ioerr.eq.0) then
          if (numf.ne.nfield) ifail = 1
        else
          read(numfs(1:4),*,iostat=ioerr) numf
          if (ioerr.eq.0) then
            if (numf.ne.nfield) ifail = 1
          else
            ifail=1
          end if
        end if
      end if
      
      if (ptra.ge.0) then
      itest=0

c *** set crude boundaries around the given point to reduce search time
      call bounds (ptra,ptdec,lora,hira,lodec,hidec)
      if (debug) write(6,*)ptra,ptdec,lora,hira,lodec,hidec

c *** convert given ra,dec (hhmmt,#ddmm) to decimal values 
      call convrt (ptra,ptdec,sra,sdec)
      if (debug) write(6,*)ptra,ptdec,sra,sdec

c *** check to see if plate centre is within crude proximity of the
c *** specified point (search has to allow for lora > hira, i.e.
c *** when range lora to hira passes through 0 hours)
      if ((((hira.gt.lora).and.(ira.ge.lora).and.(ira.le.hira)).or.
     : ((hira.lt.lora).and.((ira.le.hira).or.(ira.ge.lora)))).and.
     :    (idec.ge.lodec).and.(idec.le.hidec)) then

c *** convert plate centre ra,dec (hhmmt,#ddmm) to decimal values
         call convrt (ira,idec,pra,pdec)
      if (debug) write(6,*)ira,idec,pra,pdec

c *** following code converts ra,dec to standard co-ordinates
         radiff = sra-pra
         if (radiff.lt.-12.0) radiff=radiff+24.0
         if (radiff.gt.+12.0) radiff=radiff-24.0
         radiff=pi/12.0 * radiff
         sdeca=pi/180.0 * sdec
         pdeca=pi/180.0 * pdec
         call stanco (radiff,pdeca,sdeca,xi,eta)
         xpl=xi*c
         ypl=eta*c
      if (debug) write (6,*)xpl,ypl

c *** check stand. coords of plate centre within range of specified point
         if (abs(xpl).le.177.8.and.abs(ypl).le.177.8) itest=1
      endif
      if (itest.ne.1) ifail=1
      end if

      if (ifail.eq.0) then
c        write (out2,'(a)') line
         read (line,'(A8,A7,3A5,A6,A4,2A6,A4,A8)',iostat=ioerr)
     +        s1,s2,s3,s4,s5,s6,s7,s8,s9,s10,s11
         if (ioerr.ne.0) call all_crap('reading line for output?')
         write(out2,'(X,A8,X,A7,2X,A5,2X,A5,X,A5,2X,A6,X,A4,X,A6,X,A6,2X,
     +        A4,2X,A8)') 
     +        S1,S2,S3,S4,S5,S6,S7,S8,S9,S10,S11

        icount = icount+1
      end if
      if(debug.and.nread.gt.100) goto 999
      goto 204

 998  call all_crap('error reading catalogue')
      stop
 999  write(out2,'(x,''</PRE>'')')
      write(out2,'(x,''<P><H3><B>'')')
      write(out2,'(x''Number of plates read  = ''i8'' <br>'')') nread
      write(out2,'(x''Number plates chosen  = ''i8)') icount
      write(out2,'(x,''</H3></B>'')')
      write(out2,'(x,''<P>'')')
      write(out2,'(x,''Press "BACK" to refine this search or'')')
      write(out2,'(x,''<a href="http://www.roe.ac.uk/ukstu/''
     +  ''ukscat.html" >go to a new search.</a>'')')
      write(out2,'(x,''</HTML>'')')

      end

c      subroutine all_crap(message)
c      character*(*)  message
c      write (6,'(a)') message
c      end

      subroutine bounds (ira,idec,lora,hira,lodec,hidec)

c+
c   given
c   -----
c       ira - specified point ra (hhmmt)
c       idec - specified point dec (#ddmm)
c
c   returns
c   -------
c       boundaries for search of plates :-
c          lora - hira  for right ascension
c          lodec - hidec for declination
c
c   bounds are set as follows
c            |dec| < 70     plus or minus 1 hour  ra, 4 degrees dec
c       70 < |dec| < 79     plus or minus 3 hours ra, 4 degrees dec
c       79 < |dec|        all ras, all decs from 75 degrees to pole
c
c-
      integer ira,idec,lora,hira,lodec,hidec
      integer dec,ind,jra(2)
      data jra/1000,3000/

      dec=abs(idec)

c   if dec less than 70 degrees

      if (dec.le.7000) then
         ind=1

c   if dec between 70 and 79 degrees

      elseif ((dec.gt.7000).and.(dec.le.7900)) then
         ind=2

c   if dec greater than 79 dgrees (i.e. near the poles)

      else
         ind=3
      endif

c   if not near the poles

      if (ind.ne.3) then

c   set approximate ranges

         lora=ira-jra(ind)
         hira=ira+jra(ind)
         lodec=idec-400
         hidec=idec+400

c   perform corrections to allow for ranges crossing 0 hrs. ra

         if (lora.lt.0) lora=lora+24000
         if (hira.gt.24000) hira=hira-24000

c   if near the poles

      else
         
c   set maximum range for ra

         lora=0
         hira=24000

c   set declination range according to which hemisphere

         if (idec.lt.0) then
            lodec=-9000
            hidec=-7500
         else
            lodec=7500
            hidec=9000
         endif

      endif

      end


      subroutine convrt (ira,idec,dpra,dpdec)

c+
c   given
c   -----
c       ira  - integer value of ra  (hhmmt)
c       idec - integer value of dec (#ddmm)
c
c   returns
c   -------
c       dpra  - double precision value for hours (ra)
c       dpdec - double precision value for degrees (dec)
c
c-

      integer ira,idec,hours,adegs,degs,mins,tmins
      double precision dpra,dpdec,fhours,fdegs,ddec

c   ra conversion

      hours=ira/1000
      tmins=ira-hours*1000
      fhours=dble(tmins)/6d02
      dpra=dble(hours)+fhours

c   dec conversion

      ddec=dble(idec)
      adegs=abs(idec)
      degs=adegs/100
      mins=adegs-degs*100
      fdegs=dble(mins)/6d01
      dpdec=dsign((fdegs+dble(degs)),ddec)

      end

      subroutine stanco(radiff,pldec,stdec,xi,eta)
c+
c           computes standard coordinates for
c           a given 'star' position and plate
c           centre.   requires (in radians)
c           radiff, the difference in ra between
c           plate centre and star (star minus
c           plate), and pldec and stdec, the
c           declinations of plate centre and
c           'star' respectively.   returns
c           standard coordinates xi and eta
c           in radians.
c-
      double precision radiff,pldec,stdec,xi,eta,cds,sds,cdp,sdp,cr,sr,z
c
      cds = cos(stdec)
      sds = sin(stdec)
      cdp = cos(pldec)
      sdp = sin(pldec)
      cr  = cos(radiff)
      sr  = sin(radiff)
c
      z   = cr*cds*cdp + sds*sdp
      xi  = sr*cds/z
      eta = (sds*cdp - cr*cds*sdp)/z
c
      return
      end


      subroutine all_crap( error_message )

      character*(*) error_message

      write(6,*) '<HTML>'
cccc      write(6,*) '<BODY BACKGROUND="http://www.aao.gov.au/local/www/surveys/apm/forms/yellow.gif">'

      write(6,*) '<CENTER>'
c      write(6,*) '<IMG SRC="http://www.aao.gov.au/local/www/surveys/apm/forms/transparentfile.gif" '//
c     &     "></a>" 
      write(6,*) '<h2> Sorry, there seems to be something wrong:</h2>'
      
      write(6,*) '<p>'

      write(6,*) '<h2> ',error_message,' </h2>'

      write(6,*) '<p>'

      write(6,*) '<h2> Please try again.</h2>'
      write(6,*) '</CENTER>'

      write(6,*) '</html>'

      stop

      end

