Recipe: BandOntoSphere1.0

created by jeff101

Profile

• Name: BandOntoSphere1.0
• ID: 102961
• Shared with: Public
• Parent: None
• Children: None
• Created on: December 12, 2018 at 05:50 AM UTC
• Updated on: December 12, 2018 at 05:50 AM UTC
• Description:

  *1.0a3.txt 12/11/18 1132pm code

Code

-- This recipe takes the end point of the most 
-- recent band to space (bhub) and uses it as the 
-- center of a sphere. It then adds nb bands to the 
-- center of the sphere from alpha-carbons in the 
-- protein at segments nlo to nhi step nseg. Then, 
-- if keephub==0, it removes the original band 
-- (bhub) that pointed to the center of the sphere.
--
-- If it can't find any bands to space when it 
-- starts, this recipe will use a zero-length band 
-- from alpha-carbon (nlo+nhi)/2 as bhub to mark the 
-- center of the sphere.
--
-- As it stands, this recipe's bands will make the 
-- protein move onto the surface of the sphere.
-- If 1-way bands were used instead, 
-- ones that keep band length < goallength
-- would make this recipe be BandIntoSphere,
-- and ones that keep band length > goallength
-- would make this recipe be BandOutOfSphere.
-- 
-- Begun 12/11/18 415am by jeff101.
-- Last updated 12/11/18 1132pm.

function main()
rname='BandOntoSphere1.0' -- recipe name
nres=structure.GetCount() -- # of residues or segments
nbands=band.GetCount() -- # of bands already in use
oldscore=current.GetScore() -- present score of the structure
bhub=0 -- band # for hub band

-- below tries to find a good hub band bhub
-- from the nbands bands already in use.
for bnum=nbands,1,-1 do -- bnum counts backwards from nbands to 1
   if bhub==0 then
      tryb=band.AddToBandEndpoint(1,bnum) -- try band CA1 to end of bnum
      if tryb then -- tryb worked
         bhub=bnum
         band.Delete(tryb)
         print('Will center sphere on end point of hub band bhub='..bhub..'.')
      end -- if tryb
   end -- if bhub
end -- for bnum

-- below are default values for inputs
nlo=1
nhi=nres
nseg=7
blen=18
bstr=1
keephub=1

-- below reads the inputs
dlog=dialog.CreateDialog(rname)
dlog.nlo=dialog.AddSlider("nlo",nlo,1,nres,0)
dlog.nlolbl=dialog.AddLabel("nlo: lowest segment # to band to bhub")
dlog.nhi=dialog.AddSlider("nhi",nhi,1,nres,0)
dlog.nhilbl=dialog.AddLabel("nhi: highest segment # to band to bhub")
dlog.nseg=dialog.AddSlider("nseg",nseg,1,nres,0)
dlog.nseglbl=dialog.AddLabel("nseg: # of segments between bands")
dlog.blen=dialog.AddSlider("blen",blen,0,100,0)
dlog.blenlbl=dialog.AddLabel("blen: band goal length or sphere radius")
dlog.bstr=dialog.AddSlider("bstr",bstr,0,10,1)
dlog.bstrlbla=dialog.AddLabel("bstr: band strength")
dlog.bstrlblb=dialog.AddLabel("strong bands give a thin shell")
dlog.bstrlblc=dialog.AddLabel("weak bands give a thick shell")
dlog.keephub=dialog.AddSlider("keephub",keephub,0,1,0)
dlog.keephublbl=dialog.AddLabel("keephub: 1 keeps hub band bhub, 0 removes it")
dlog.run=dialog.AddButton("Run",1)
dlog.quit=dialog.AddButton("Quit",0)
dlogcode=dialog.Show(dlog)

if dlogcode>0 then
   -- below stores inputs
   nlo=dlog.nlo.value
   nhi=dlog.nhi.value
   nseg=dlog.nseg.value
   blen=dlog.blen.value
   bstr=dlog.bstr.value
   keephub=dlog.keephub.value

   -- below checks inputs
   if nlo>nhi then
      nlo,nhi=nhi,nlo -- swap em
   end -- if nlo

   -- below prints the inputs
   print(string.format('Got the inputs nlo=%d nhi=%d nseg=%d', 
      nlo,nhi,nseg))
   print(string.format('blen=%d bstr=%.1f keephub=%d.',
      blen,bstr,keephub))

   if bhub==0 then -- still need to make a hub band #'d bhub
      nmid=round(0.5*(nlo+nhi)) -- pick middle segment of nlo to nhi

      -- below picks 2 other segments nx,ny to set up 
      -- band to space coordinate system with nmid as its origin.
      isok=0
      nx,ny=1,nres
      incx=1 -- 1 means add to nx, 0 means subtract from ny
      while(isok==0) do -- do until isok is not zero
         if nmid~=nx and nx~=ny and ny~=nmid then
            d1=getdist(nmid,nx)
            d2=getdist(nx,ny)
            d3=getdist(ny,nmid)

            -- below gets the order 0<=d1<=d2<=d3
            if d1>d2 then
               d1,d2=d2,d1 -- swap em
            end -- if d1
            if d1>d3 then
               d1,d3=d3,d1 -- swap em
            end -- if d1
            if d2>d3 then
               d2,d3=d3,d2 -- swap em
            end -- if d2
            if d1>=0 then
               if math.abs(d1+d2-d3)>0.01 then
                  -- nmid,nx,ny are not collinear, a good thing
                  isok=1 -- this will quit the loop!
               end -- if math
            end -- if d1
         end -- if nmid
         if isok==0 then -- try to pick a new set of nx,ny
            if incx==1 then -- add to nx
               nx=nx+1
               if nx>nres then
                  print('Error: nx='..nx..' is > nres='..nres..'.')
                  return
               end -- if nx
            else -- subtract from ny
               ny=ny-1
               if ny<1 then
                  print('Error: ny='..ny..' is > nres='..nres..'.')
                  return
               end -- if ny
            end -- if incx
            incx=1-incx -- have incx swap between 0 and 1
         end -- if isok
      end -- while isok 

      -- below tries to make the hub band
      tryb=band.Add(nmid,nx,ny,0.001,0,0)
      -- Above uses alpha-carbons (CA) on residues nmid,nx,ny 
      -- to set up xyz coordinates where nmid's CA is the origin, 
      -- nx's CA is along the +x axis, and ny's CA sets up the xy plane 
      -- with ny's CA having y>0. The band to space will extend a distance
      -- 0.001 A from nmid's CA along the +z axis.
      --
      -- 12/11/18 Can't have segment for origin & +x axis be the same
      -- even if using different atom #'s. Using CA,C,N on nmid would
      -- have been easier if it were allowed.
      --
      if tryb then -- tryb worked
         band.SetStrength(tryb,1)
         band.SetGoalLength(tryb,0)
         band.Disable(tryb)
         bhub=tryb
         print('Will center sphere on end point of hub band bhub='..bhub)
         print('(very close to alpha-carbon '..nmid..' at first).')
      else
         print('Error making zero-length band from alpha-carbon '..nmid..'.')
      end -- if tryb
   end -- if bhub

   -- below uses & resets some inputs
   -- In below,
   -- if nlo=1 & nhi=25 and have nseg=8,
   -- want nb=4 bands from segments 1 9 17 25 to center.
   -- Also, if nlo=1 & nhi=37 and have nseg=6,
   -- want nb=7 bands from segments 1 7 13 19 25 31 37 to center.
   nb=1+math.ceil((nhi-nlo)/nseg)
   nseg=(nhi-nlo)/(nb-1)

   -- below lists inputs
   mssg=string.format('Adding %d bands spaced %d-%d segments apart with',
      nb,math.floor(nseg),math.ceil(nseg))
   mssg=string.format('%s\ngoal length %d A and strength %.1f',mssg,blen,bstr)
   mssg=(mssg..' to the end point of\nhub band bhub='..bhub)
   mssg=(mssg..' from the alpha-carbons of segments\n')

   numgot=0 -- total times tryb has worked so far
   -- below makes nb new bands
   for num=1,nb do 
      roundnum=round(nlo+(num-1)*nseg)
      tryb=band.AddToBandEndpoint(roundnum,bhub)
      if tryb then -- tryb worked
         numgot=numgot+1
         if numgot>1 then        
            mssg=(mssg..' '..roundnum)
         else 
            mssg=(mssg..roundnum)
         end -- if numgot
         band.Enable(tryb)
         band.SetGoalLength(tryb,blen)
         band.SetStrength(tryb,bstr)
      end -- if tryb
   end -- for num
   print(mssg..'.')
   -- ideally, numgot=nb now
   -- and segments nlo & nhi are both listed in mssg

   -- below removes hub band bhub if it exists and keephub==0
   mssg=('keephub='..keephub)
   if keephub==0 then
      if bhub~=0 then
         print(mssg..', so removing hub band bhub='..bhub..'.')
         band.Delete(bhub) -- remove the original band 
         -- to the center of the sphere. all band indices above 
         -- bhub now shift down by 1 to compensate.
      else 
         print(mssg..', but there is no hub band to remove.')
      end -- if bhub
   else
      print(mssg..', so keeping any hub band that exists.')
   end -- if keephub

   -- below gives final message
   print(string.format('Began with nres=%d nbands=%d and score=%.3f',
      nres,nbands,trunc3to0(oldscore)))
   print(string.format('and ending with nbands=%d and score=%.3f.',
      band.GetCount(),trunc3to0(current.GetScore())))
end -- if dlogcode
end -- main ()

-- below finds the distance between the 
-- alpha-carbons on segments seg1 & seg2.
function getdist(seg1,seg2)
local vres= -1 -- the distance
local bnum -- # of a dummy band used to measure distances
local segdiff=math.abs(seg1-seg2) -- |seg1-seg2|
if segdiff>1 then
   vres=structure.GetDistance(seg1,seg2)
elseif segdiff==1 then
   bnum=band.AddBetweenSegments(seg1,seg2) 
   -- bnum is a dummy address
   if bnum then -- band bnum was created fine
      vres=band.GetLength(bnum)   
      band.Delete(bnum)
   end -- if bnum
end 
if vres<0 then
   print(string.format('Error: distance between alpha-carbons %d and %d is %.1f.',
      seg1,seg2,vres))
end -- if vres
return vres
end -- getdist()

-- below rounds the number val
-- to the nearest integer vres
function round(val)
local vres=math.floor(val)
if val-vres>=0.5 then
   vres=vres+1
end -- if vres
return vres
end -- round()

-- below is how Foldit truncates scores to the 0.001 position:
-- for positive scores, the score drops, and
-- for negative scores, the score rises, 
-- either way truncating the score towards 0.
function trunc3to0(val)
return val-(val%0.001)
end -- trunc3to0()

main()