implemented printPDBnotfound option

Author: lmiq

input/inputfile.inp

@@ -16,7 +16,6 @@ printlinks no
 linkdir ./links
 
 # What cross-links to compute (all/reactive/observed)
-
 compute observed
 #compute reactive
 #compute all
@@ -99,14 +98,15 @@ end experiment DSS
 
 search_limit relative 1.5  # Search for topological distances of up to search_range*dmaxlink
 
-printallfound no  # Prints, or not, the structure of the links that were found but violate the distances
+printallfound no    # Prints, or not, the structure of the links that were found but violate the distances
+printPDBnotfound no # Prints, or not, PDB files for the links which were not found at all (prints a straight line)
 
 pgood 0.70 # Probability of observing a link which is within link reach
 pbad  0.01 # Probability of observing a link which is NOT within linker reach
 
 vdwrad 3.         # Radius to be used for excluded atom volume computation
 print 0           # Extensive or concise printing
-printnotfound yes # Print data for pairs for which links were not found
+printnotfound yes # Print data for pairs for links which were not found
 quitgood no       # Quit whenever a valid link is found (stop searching for shortest one)
 dbond 1.5         # Bond distance of link "atoms"
 kbond 10.         # Energy constant for link bonds

src/inputoptions.f90

@@ -29,6 +29,7 @@ module inputoptions
   logical :: printlinks
   logical :: printallfound
   logical :: printnotfound
+  logical :: printPDBnotfound
   logical :: observedscores
   logical :: mimicchain
   logical :: printaccessible

src/link_to_pdb.f90

@@ -10,8 +10,8 @@ subroutine link_to_pdb(link,nlinkatoms,linkdir,pdbfile,natoms,atom,n,x)
   implicit none
 
   integer :: ioerr
-  integer :: j, ix, iy, iz
-  integer :: natoms, nlinkatoms
+  integer :: j, ix, iy, iz, index
+  integer :: natoms, nlinkatoms, nsteps
 
   type(specific_link) :: link
   type(pdbatom) :: writeatom, atom(natoms)
@@ -20,7 +20,7 @@ subroutine link_to_pdb(link,nlinkatoms,linkdir,pdbfile,natoms,atom,n,x)
   character(len=4) :: char1, char2
 
   integer :: n
-  double precision :: x(n)
+  double precision :: x(n), stepx, stepy, stepz
   double precision :: overlap, stretch
   character(len=13) :: statuschar
 
@@ -41,33 +41,57 @@ subroutine link_to_pdb(link,nlinkatoms,linkdir,pdbfile,natoms,atom,n,x)
     stop
   end if
   write(11,"( 'REMARK EUCLIDEAN DISTANCE: ', (tr2,f12.5) )") link%euclidean
-  write(11,"( 'REMARK TOPOLOGICAL DISTANCE: ', (tr2,f12.5) )") link%topodist
-  write(11,"( 'REMARK OVERLAP and STRETCH: ', 2(tr2,f12.5) )") overlap(n,x), stretch(n,x)
   write(11,"( 'REMARK LINK STATUS: ', a13 )") statuschar(link%status)
+  if ( link%found ) then
+    write(11,"( 'REMARK TOPOLOGICAL DISTANCE: ', (tr2,f12.5) )") link%topodist
+    write(11,"( 'REMARK OVERLAP and STRETCH: ', 2(tr2,f12.5) )") overlap(n,x), stretch(n,x)
+  end if
+  index = 1
   writeatom = atom(link%atom1%index)
   writeatom%residue%name = "LINK"
   writeatom%residue%chain = "A"
   writeatom%residue%index = 1
-  writeatom%index = 1
+  writeatom%index = index
   write(11,"(a)") trim(print_pdbhetatm(writeatom))
-  do j = 1, nlinkatoms
-    ix = (j-1)*3 + 1
-    iy = ix + 1
-    iz = ix + 2
-    writeatom%index = j + 1
-    writeatom%name = "O"
-    writeatom%residue%name = "LINK"
-    writeatom%residue%chain = "A"
-    writeatom%residue%index = 1
-    writeatom%x = x(ix)
-    writeatom%y = x(iy)
-    writeatom%z = x(iz)
-    write(11,"(a)") trim(print_pdbhetatm(writeatom))
-  end do
+  if ( link%found ) then
+    do j = 1, nlinkatoms
+      index = index + 1
+      ix = (j-1)*3 + 1
+      iy = ix + 1
+      iz = ix + 2
+      writeatom%index = index
+      writeatom%name = "N"
+      writeatom%residue%name = "LINK"
+      writeatom%residue%chain = "A"
+      writeatom%residue%index = 1
+      writeatom%x = x(ix)
+      writeatom%y = x(iy)
+      writeatom%z = x(iz)
+      write(11,"(a)") trim(print_pdbhetatm(writeatom))
+    end do
+  else
+    nsteps = int(link%euclidean)-2
+    stepx = (link%atom2%x - link%atom1%x)/nsteps
+    stepy = (link%atom2%y - link%atom1%y)/nsteps
+    stepz = (link%atom2%z - link%atom1%z)/nsteps
+    do j = 1, nsteps
+      index = index + 1
+      writeatom%index = index
+      writeatom%name = "O"
+      writeatom%residue%name = "LINK"
+      writeatom%residue%chain = "A"
+      writeatom%residue%index = 1
+      writeatom%x = link%atom1%x + j*stepx
+      writeatom%y = link%atom1%y + j*stepy
+      writeatom%z = link%atom1%z + j*stepz
+      write(11,"(a)") trim(print_pdbhetatm(writeatom))
+    end do
+  end if
+  index = index + 1
   writeatom = atom(link%atom2%index)
   writeatom%residue%name = "LINK"
   writeatom%residue%chain = "A"
-  writeatom%residue%index = 1
+  writeatom%residue%index = index
   writeatom%index = nlinkatoms + 2
   write(11,"(a)") trim(print_pdbhetatm(writeatom))
   close(11)

src/printdata.f90

@@ -19,7 +19,7 @@ subroutine printdata(print,link)
 
   if ( print == -1 ) then
     write(str,dashes) ; call writelog(str)
-    write(str,"( '        RESIDUE1   ATOM1 RESIDUE2   ATOM2 EUCLDIST  TOPODIST OBSERVED&
+    write(str,"( '>     RESIDUE1   ATOM1 RESIDUE2   ATOM2 EUCLDIST  TOPODIST OBSERVED&
                  &    DMIN      DMAX        RESULT  OBSRES REACRES RA AA')") ; call writelog(str)
     write(str,dashes) ; call writelog(str)
     return

src/topolink.f90

@@ -80,6 +80,7 @@ program topolink
   printlinks = .false.
   printnotfound = .false.
   printallfound = .false.
+  printPDBnotfound = .false.
   compute = 2
   pgood = -1.d0
   pbad = -1.d0
@@ -227,6 +228,9 @@ program topolink
       case ("printnotfound")
         if ( keyvalue(record,1) == 'yes' ) printnotfound = .true.
         if ( keyvalue(record,1) == 'no' ) printnotfound = .false.
+      case ("printPDBnotfound")
+        if ( keyvalue(record,1) == 'yes' ) printPDBnotfound = .true.
+        if ( keyvalue(record,1) == 'no' ) printPDBnotfound = .false.
       case ("printaccessible")
         if ( keyvalue(record,1) == 'yes' ) printaccessible = .true.
         if ( keyvalue(record,1) == 'no' ) printaccessible = .false.
@@ -1307,6 +1311,11 @@ program topolink
           call linkconsistency(link(i),nexp,experiment)
           call printdata(print,link(i))
         end if
+        if ( printlinks ) then
+          if ( printPDBnotfound ) then
+            call link_to_pdb(link(i),nlinkatoms,linkdir,pdbfile,natoms,atom,n,xbest)
+          end if
+        end if
         cycle allpairs
       end if
 
@@ -1319,6 +1328,11 @@ program topolink
           call linkconsistency(link(i),nexp,experiment)
           call printdata(print,link(i))
         end if
+        if ( printlinks ) then
+          if ( printPDBnotfound ) then
+            call link_to_pdb(link(i),nlinkatoms,linkdir,pdbfile,natoms,atom,n,xbest)
+          end if
+        end if
         cycle allpairs
       end if
 
@@ -1456,13 +1470,20 @@ program topolink
         if ( printallfound .or. &
              ( link(i)%status == 0 .or. link(i)%status == 5 ) ) then
           call link_to_pdb(link(i),nlinkatoms,linkdir,pdbfile,natoms,atom,n,xbest)
+        else if( printPDBnotfound ) then
+          call link_to_pdb(link(i),nlinkatoms,linkdir,pdbfile,natoms,atom,n,xbest)
         end if
       end if
 
     else
 
       ! If a topological distance was not found for this pair
       if ( printnotfound ) call printdata(print,link(i))
+      if ( printlinks ) then
+        if ( printPDBnotfound ) then
+          call link_to_pdb(link(i),nlinkatoms,linkdir,pdbfile,natoms,atom,n,xbest)
+        end if
+      end if
 
     end if