	options /check=bounds
C
	program search
C
C  Search for amino acid sequence in PDB
C
C  Input: PDBSEQ, file created by MAKEDB.
C	  Looks for original Brookhaven files in directory PDBDIR:
C
C  Dan Bloch
C  Molecular Biology Dept.   MB13
C  Research Institute of Scripps Clinic
C  La Jolla, CA  92037
C
	implicit none
	character*4 protein,oldprot
	character*1 chain
	character*1000 sequence
	character*100 subseq
	integer len,sublen,ires,locate,ibase,iloc
	logical oneshot /.false./

	open(1,status='old',name='pdbseq',form='unformatted',shared,readonly)
	call lib$get_foreign(subseq,,sublen)
	if (sublen .ne. 0) then
	   oneshot = .true.
	   goto 100
	endif

1	write (*,10)
10	format (/' Sequence:  ',$)
	read (*,20,end=999) sublen,subseq
20	format (q,a)
	if (sublen.eq.0.or.subseq.eq.' ') goto 999
	call str$upcase(subseq,subseq)
	rewind 1

100	do while (.true.)
	   read (1,end=888) protein,chain,len,sequence(1:len)
	   ibase=0
	   do while (ibase+sublen.le.len)
	      iloc = locate(sequence(ibase+1:len),subseq(1:sublen))
	      if (iloc.eq.0) then
		 ibase=len+1	! so we exit do loop
	      else
	         if (protein .ne. oldprot) then
		    call identify(protein)
	            oldprot=protein
		 endif
		 ires=ibase+iloc
	         if (chain .eq. ' ') then
	            write (*,110) sequence(ires:ires+sublen-1),ires
110		    format(4x,a,' starting at residue',i4,:,' in chain ',a1)
	         else
	            write (*,110) sequence(ires:ires+sublen-1),ires,chain
	         endif
		 ibase=ires+1
	      endif
	   enddo
	enddo

888	if (.not. oneshot) goto 1
999	end
C-----------------------------------------------------------------------------
	subroutine identify(protein)
C
C  Give PDB identifier, type out COMPND and SOURCE lines from file.
C
	implicit none
	character*(*) protein
	character*80 line
	character*6 type
	equivalence (line,type)
C
	open(11,status='old',shared,readonly,err=66,
	1	name='PDB'//protein,defaultfile='PDBDIR:.ENT')

	write (*,*)
	type = 'HEADER'
	do while (type.eq.'HEADER')
	   read(11,110)line
110	   format (a80)
	enddo

	do while (type.eq.'COMPND'.or.type.eq.'SOURCE')
	   write (*,120) line
120	   format (1x,a80)
	   read(11,110)line
	enddo
	close (11)
	return

66	write (*,70) protein
70	format (/' Protein not found: ',a)
	return
	end
C------------------------------------------------------------------------------
	function locate(string1,string2)
C
C  Amino acid string search with B,Z,X wildcards.  Syntax is like Fortran INDEX.
C
	character*(*) string1,string2
	character*1 char1,char2

	len1=len(string1)
	len2=len(string2)
	do 20 i=0,len1-len2
	   do 10 j=1,len2
	      char1=string1(i+j:i+j)
	      char2=string2(j:j)
	      if (char2.eq.'X'.or.char2.eq.'*') goto 10
	      if (char2.eq.'B'.and.(char1.eq.'N'.or.char1.eq.'D')) goto 10
	      if (char2.eq.'Z'.and.(char1.eq.'Q'.or.char1.eq.'E')) goto 10
	      if (char2.eq.char1) goto 10
	      goto 20
10	   continue
	   locate=i+1	! success
	   return
20	continue

33	locate=0	! failure
	return
	end
