Recipe: V2 multibander v3.13

created by jeff101

Profile

• Name: V2 multibander v3.13
• ID: 44555
• Shared with: Public
• Parent: None
• Children: None
• Created on: December 15, 2012 at 15:50 PM UTC
• Updated on: December 15, 2012 at 15:50 PM UTC
• Description:

  Adapted from V2 multibander v3.12 11/10/12.

Code

-- MB
-- multibander v3.12 last updated 11/4/12 1210am midnight
-- multibander v3.13 begun 1038pm 11/10/12 and last updated 1235am midnight 11/11/12

-- calculate the number of ss changes for LS Quake
cnt=0 lastss="?" for seg=1,structure.GetCount() do if structure.GetSecondaryStructure(seg)~=lastss then cnt=cnt+1 lastss=structure.GetSecondaryStructure(seg) end end

print(string.format("Running V2 multibander v3.13 @ %s",os.date()))

ped={}
ped.use_gui=true
ped.bander_order={}
ped.bands_at_start=band.GetCount()
print('Starting with '..ped.bands_at_start..' bands.')

--ped.bander_order="QS,C,QR,VK,LSQ" keep in for reference
--  "QS"  -> "Quake Slice"
--  "C"   -> "Compressor"
--  "QR   -> "QuakeR"
--  "VK   -> "Voids Killer"
--  "LSQ" -> "LS Quake"

ped.bander_order="QS,C,QR,VK,LSQ"
ped.run_quake_slice=true
ped.run_compressor=true
ped.run_quake=true
ped.run_voids_killer=true
ped.run_ls_quake=true

ped.quake_slice_maxns=5
ped.quake_slice_un=true
ped.quake_slice_deux=true
ped.compressor_loops=50
ped.compressor_modulator=true --true --changing compression/decompression each loop
ped.compressor_allLoop=false --false --work in all-loop mode. sometimes work better than structure mode :)
ped.quake_loops=50
ped.voids_killer_loops=50
ped.ls_quake_loops=3*cnt      -- default to 3 times the number of ss changes
ped.total_loops=10
ped.relax_protein_option=1
--                          0 - no relaxation of protein
--                          1 - wiggle down
--                          2 - shake down
--                          3 - rebuild 10 worst segments
--                          4 - pull hydrophobes and push hydrophilics
--                          5 - rebuild the whole protein
ped.minimum_loop_gain=10    -- minimum loop gain that if not met will trigger relaxing the protein (set to 9999 to always relax protein)
ped.minimum_gain1=1         -- minimum gain that bander needs to meet to avoid being temporarily disabled
ped.minimum_gain2=1         -- minimum gain that needs to be made before bander will be re-enabled
ped.relax_from_best=true    -- restore best score before relaxing protein
ped.relaxation_minimum_amount=100 -- minimum amount of points needed to be lost for relax_protein_option=3
ped.enable_mutate=true      -- mutate after every bander if mutation puzzle
ped.bcis={}
--ped.bcis={0.1,0.05,0.01}
ped.bcis={0.005,0.01,0.015}
ped.bci_index=1
ped.gains={}
ped.options={}
ped.options={"bander_order","run_quake_slice","run_compressor","run_quake","run_voids_killer","run_ls_quake","separator",
             "quake_slice_maxns","quake_slice_un","quake_slice_deux","separator","compressor_loops","compressor_modulator",
             "compressor_allLoop","separator","quake_loops","separator","voids_killer_loops","separator","ls_quake_loops","separator","total_loops","relax_protein_option",
             "minimum_loop_gain","minimum_gain1","minimum_gain2","relax_from_best","relaxation_minimum_amount","enable_mutate"}

count=0
for seg=1,structure.GetCount() do
  if structure.IsMutable(seg) then count=count+1 end
end
if count==0 then ped.enable_mutate=false end

-- =========  change history =========================================
-- multibander v3.13
-- make output scores match Foldit GUI's scores using trunc3()
-- list rms diff for user-supplied bands in listdiff()
-- list which AA's get mutated

-- multibander v3.12 
-- fixed Void Killer bug noticed in puzzle 644 for non-amino-acid segment 31

-- multibander v3.11
-- adapted from v3.10 to not delete bands present at the start
 
-- multibander v3.10
-- added gui interface

-- multbander v3.09
-- added relax_protein_option 5

-- multbander v3.08
-- fixed problem with mutating code

-- multbander v3.07
-- added recording information in Note

-- multbander v3.06
-- added mutating for mutation puzzles

-- multbander v3.05
-- fixed function KillVoids (se=math.min(se,structure.GetCount()))

-- multibander v3.04
-- added push_and_pull routine option to ped.relax_protein_option

-- multibander v3.03
-- attempting to include updates from Compressor v2.2

-- multibander v3.02
-- added bestScore=Score() in Voids Killer, needed for 1st run after relaxing protein
-- initialized ped.disable["LS Quake"] to false

-- multibander v3.01
-- added programmable ordering of banders (ped.bander_order)
-- added LS Quake

-- multibander v2.04
-- added bestScore=Score() in quake slice, needed for 1st run after relaxing protein

-- multibander v2.03
-- added new option 3 for relaxing protein (don't use it though)
-- print options at start
-- replaced quake slice v1.2 with quake slice v1.3

-- multibander v2.02
-- fixed make_a_decision for disabling bander

-- multibander v2.01
-- resolved function conflict by adding _MWC2 function

-- multibander v2.0
-- put all functions at top
-- removed duplicate functions (resolved differences in functions)
-- added option to conditionally relax protein at the end of the loop
-- added conditions to disable and re-enable banders

-- multibander v1.0
-- four banders run one after another for 50 loops, options at top
-- first bander - Quake slice v1.3
-- second bander - Rav3n_pl Compressor v2.0
-- third bander - QuakeR  v3.0
-- fourth bander - Rav3n_pl Voids Killer v0.4
-- fifth bander - LS Quake
-- =============================================================================
--[[
Quake slice - a Quake by slice by Marie Suhard
based on Rav3n Quake R
code cleaned by rav3n_pl
]]--

--[[
    Rav3n_pl Compressor v2
    trying to compress/decompress protein accepting loss of points between pulls
]]--

--[[
QuakeR - a randomized Quake
Last update 2011-04-18 rav3n_pl

Based on "Quake" by Grom
]]--

--[[
Voids Killer by rav3n_pl
script is searching for possible banding across voids.
Because we not have any voids detection in Lua, I do some math to found them
Not always perfect, but mostly works :)
]]--

--[[
LS Quake
Base on Local Quake 1.0 by spvinvent  with helep by thom
Bands radiate from a structure (helix, loop or sheet) to spatially close neighbours
one round with the strcuture frozen one round not
Foregoes Fuses and generally tries to optimize for speed
]]--
-- =============================================================================

-- below rounds val down to nearest 0.001
function trunc3(val)
return val-(val % 0.001)
end -- trunc3()

-- below rounds to 2 decimal places
function round2(numi)
local numo=round(numi*100)/100
return numo
end -- round2()

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

function listdiff()
-- new to v3.13. adapted from adjustbands2d3.txt
-- lists how well present conformation obeys user-supplied bands
local nbands, lens, glens, strs, diffs, ison, tots, avgdiff, rmsdiff, btype, mindiff, maxdiff
local maxstr, minstr, num, str, diff, diff2, indx, btype, typenams
nbands=ped.bands_at_start
lens={}
glens={}
strs={}
diffs={}
ison={}
tots={}
avgdiff={}
rmsdiff={}
for btype=0,2 do
   avgdiff[btype]=0
   rmsdiff[btype]=0
   tots[btype]=0
end -- for btype
mindiff={}
maxdiff={}
maxstr={}
minstr={}
for num=1,nbands do
   ison[num]=band.IsEnabled(num)
   lens[num]=band.GetLength(num)
   glens[num]=band.GetGoalLength(num)
   str=band.GetStrength(num)
   strs[num]=str
   diff=math.abs(lens[num]-glens[num])
   diffs[num]=diff
   diff2=diff*diff
   if ison[num]==true then 
      btypes={1,2} -- for enabled bands
   else 
      btypes={0,2} -- for disabled bands
   end -- if ison
   for indx=1,2 do
      btype=btypes[indx]
      avgdiff[btype]=avgdiff[btype]+diff
      rmsdiff[btype]=rmsdiff[btype]+diff2
      tots[btype]=tots[btype]+1
      if tots[btype]==1 then
         mindiff[btype]=diff
         maxdiff[btype]=diff
         maxstr[btype]=str
         minstr[btype]=str
      else
         if diff>maxdiff[btype] then
            maxdiff[btype]=diff
         elseif diff<mindiff[btype] then
            mindiff[btype]=diff
         end
         if str>maxstr[btype] then
            maxstr[btype]=str
         elseif str<minstr[btype] then
            minstr[btype]=str
         end
      end -- if tots
   end -- for indx
end -- for num
for btype=0,2 do
   avgdiff[btype]=avgdiff[btype]/tots[btype]
   rmsdiff[btype]=math.sqrt(rmsdiff[btype]/tots[btype])
end
print('For '..nbands..' user-supplied bands and score '..trunc3(current.GetScore())..' we have:')
typenams={}
typenams[0]='disabled'
typenams[1]=' enabled'
typenams[2]='any-type'
for btype=0,2 do
   if tots[btype]>0 then
      print('   '..tots[btype]..' '..typenams[btype]..' bands with diffs='..round2(mindiff[btype])..'-'..round2(maxdiff[btype])..', avgdiff='..round2(avgdiff[btype])..', rmsdiff='..round2(rmsdiff[btype])..', and strengths='..round2(minstr[btype])..'-'..round2(maxstr[btype])..'.')
   end -- if tots
end -- for btype
end -- listdiff()

function getseq()
-- new to v3.13. adapted from eke3d2.txt.
-- stores present amino acid sequence in variable seq
local seq={}
local num
for num=1,structure.GetCount() do
   seq[num]=structure.GetAminoAcid(num)
end -- for num
return seq
end -- getseq()

function listseqdiffs(oldseq,newseq)
-- new to v3.13. adapted from eke3d2.txt
-- lists differences between amino acid sequences oldseq and newseq
local num
local tot=0
local str=' '
for num=1,structure.GetCount() do
   if oldseq[num]~=newseq[num] then
      tot=tot+1
      str=(str..oldseq[num]..num..'->'..newseq[num]..num..' ')
   end -- if oldseq
end -- for num
if tot>0 then
   print('   Mutated residues'..str..' ('..tot..' in all).')
else
   print('   No mutations succeeded.')
end -- if tot
end -- listseqdiffs()

function MyBandDeleteAll()
-- new to v3.11 
-- use this instead of band.DeleteAll() throughout the code
local nbands=band.GetCount()
local oldbands=ped.bands_at_start
if nbands>oldbands then
if oldbands==0 then
   print('Deleting all '..nbands..' bands.')
   band.DeleteAll() -- delete all the bands at once
else
   local bnum
   local str=('Deleting bands '..(oldbands+1)..'-'..nbands)
   for bnum=nbands,(oldbands+1),-1 do
      band.Delete(bnum) -- delete bands one by one from nbands down to oldbands+1
   end -- for bnum
   nbands=band.GetCount()
   print(str..' and leaving '..nbands..' original bands intact.')
   if nbands~=oldbands then
      print('ERROR in MyBandDeleteAll: ending with nbands='..nbands..' different from oldbands='..oldbands..'.')
   end -- if nbands
end -- if oldbands
end -- if nbands
end -- MyBandDeleteAll()

function record_score(label)
  local pedix
  local improvement
  if ped==nil or ped.scores==nil then
    if ped==nil then ped={} end
    ped.scores={}
    ped.steps={}
    ped.loops={}
    ped.steps[0]=label
    ped.scores[0]=current.GetScore()
    ped.loops[0]=0
    print(string.format("%s: %.3f",ped.steps[0],trunc3(ped.scores[0])))
  else
    ped.steps[#ped.steps+1]=label
    ped.scores[#ped.scores+1]=current.GetScore()
    ped.loops[#ped.loops+1]=ped.loop_count
    improvement=ped.scores[#ped.scores]-ped.scores[0]
    print(string.format("%s: %.3f Total improvement: %.3f",ped.steps[#ped.steps],trunc3(ped.scores[#ped.scores]),improvement))
    print(" ")
    print("------------- Recap -------------")
    for pedix=0,#ped.scores do
      print(string.format("(%i) %s: %.3f",ped.loops[pedix],ped.steps[pedix],trunc3(ped.scores[pedix])))
    end
    print("---------------------------------")
    print(" ")
    report_time(start_clock,start_time)
    print(" ")
    print("---------------------------------")
  end
end

function make_a_decision(bander)
  local gain_since_last
  if ped.gains[bander]==nil then return end
  gain_since_last=current.GetScore()-ped.last_score[bander]
  if ped.protein_was_relaxed then
    ped.disable[bander]=false
  elseif ped.run_quake_slice and not ped.disable[bander] then
    print(bander.." gain: ",ped.gains[bander]," Gain since "..bander.." was last run: ",gain_since_last)
    if (ped.gains[bander]<ped.minimum_gain1) and (gain_since_last<ped.minimum_gain2) then
      print(bander.." failed to gain the minimum number of points to prevent from being temporarily disabled")
      ped.disable[bander]=true
    end
  elseif ped.run_quake_slice and ped.disable[bander] then
    if current.GetScore()-ped.last_score[bander] then
      print("Enough gain has been achieved to re-enable "..bander)
      ped.disable[bander]=true
    end
  end
end

function sort_by_column(a,b)
  if a[sort_column]<b[sort_column] then return true else return false end
end

function get_segment_scores_array(first,last,snap_count)
  local array1={}
  local seg
  local ix,ix2
  local attr={"clashing","packing","hiding","bonding","backbone","sidechain","reference","disulfides","other"}
  if first==nil then first=1 end
  if last==nil then last=structure.GetCount() end
  if snap_count==nil then snap_count=false else snap_count=true end
  ix2=1
  for seg=first,last do
    array1[ix2]={}
    array1[ix2][1]=seg
    array1[ix2][2]=current.GetSegmentEnergyScore(seg)
    for ix=1,#attr do
      array1[ix2][ix+2]=current.GetSegmentEnergySubscore(seg,attr[ix])
    end
    if snap_count==true then
      array1[ix2][#array1[ix2]+1]=rotamer.GetCount(seg)
    end
    ix2=ix2+1
  end
  return array1
end

function report_options()
  local i,j,option,typ,out
  for i=1,#ped.options do
    option=ped.options[i]
    if option=="separator" then
      print("-----------------------------")
    else
      typ=type(ped[option])
      if typ=="boolean" then
        if ped[option] then
          print(string.format("ped.%s=true",ped.options[i]))
        else
          print(string.format("ped.%s=false",ped.options[i]))
        end
      elseif typ=="table" then
        out=ped[option][1]
        for j=2,#ped[option] do out=out..","..ped[option][j] end
        print(string.format("ped.%s=%s",option,out))
      else
        print(string.format("ped.%s=%s",option,ped[option]))
      end
    end
  end
end

function report_time(start_clock,start_time,clock_msg,time_msg)
  local seconds,minutes,hours,days
  if clock_msg==nil then clock_msg="CPU time" end
  if time_msg==nil then time_msg="Elapsed time" end
  print(string.format("%s",os.date()))
  days,remainder=math.modf((os.clock()-start_clock)/(24*60*60))
  hours,remainder=math.modf(remainder*24)
  minutes,remainder=math.modf(remainder*60)
  seconds,remainder=math.modf(remainder*60)
  print(string.format("%s(%02id:%02ih:%02im:%02is)",clock_msg,days,hours,minutes,seconds))
  days,remainder=math.modf(os.difftime(os.time(),start_time)/(24*60*60))
  hours,remainder=math.modf(remainder*24)
  minutes,remainder=math.modf(remainder*60)
  seconds,remainder=math.modf(remainder*60)
  print(string.format("%s(%02id:%02ih:%02im:%02is)",time_msg,days,hours,minutes,seconds))
end

function push_and_pull()
  local aminoacid,lastseg,maxbands,hydrophobes,hydrophilics
  local seg,aa,seg1,seg2,dist,b,score1,score2,lastband,b,ix

  aminoacid={}
  aminoacid.hydrophobes="ivlfcma"
  aminoacid.hydrophilics="rkndeqh"

  math.randomseed(os.clock())

  lastseg=structure.GetCount()
  maxbands=math.ceil(lastseg*0.1)
  hydrophobes={}
  hydrophilics={}
  for seg=1,lastseg do
    aa=structure.GetAminoAcid(seg)
    if string.find(aminoacid.hydrophobes,aa)~=nil then hydrophobes[#hydrophobes+1]=seg end
  end
  for seg=1,lastseg do
    aa=structure.GetAminoAcid(seg)
    if string.find(aminoacid.hydrophilics,aa)~=nil then hydrophilics[#hydrophilics+1]=seg end
  end

  MyBandDeleteAll()

  lastband=band.GetCount()
  for ix=1,maxbands do
    seg1=hydrophobes[math.random(#hydrophobes)]
    seg2=hydrophobes[math.random(#hydrophobes)]
    dist=structure.GetDistance(seg1,seg2)
    band.AddBetweenSegments(seg1,seg2)
    -- below checks if above line added a band successfully
    b=band.GetCount()
    if b>lastband then
      band.SetGoalLength(b,math.max(0,dist-5))
      band.SetStrength(b,0.5)
      lastband=b
    end
  end

  lastband=band.GetCount()
  for ix=1,maxbands do
    seg1=hydrophilics[math.random(#hydrophilics)]
    seg2=hydrophilics[math.random(#hydrophilics)]
    dist=structure.GetDistance(seg1,seg2)
    band.AddBetweenSegments(seg1,seg2)
    -- below checks if above line added a band successfully
    b=band.GetCount()
    if b>lastband then
      band.SetGoalLength(b,dist+5)
      band.SetStrength(b,0.5)
      lastband=b
    end
  end

  behavior.SetClashImportance(0.3+0.3*math.random())
  structure.WiggleAll(1)
  behavior.SetClashImportance(1.0)
  repeat
    score1=current.GetEnergyScore()
    structure.ShakeSidechainsAll(1)
    score2=current.GetEnergyScore()
  until score2<score1+50
  repeat
    score1=current.GetEnergyScore()
    structure.WiggleAll(2,false,true)
    score2=current.GetEnergyScore()
  until score2<score1+50
  repeat
    score1=current.GetEnergyScore()
    structure.WiggleAll(2,true,true)
    score2=current.GetEnergyScore()
  until score2<score1+50

  MyBandDeleteAll()

end

function multibander_gui()
  local rc,order
  menu=dialog.CreateDialog("multibander")
  menu.quakeslice=dialog.AddCheckbox("(QS) Quake Slice",ped.run_quake_slice)
  menu.quake_slice_maxns=dialog.AddTextbox("QS maxns",ped.quake_slice_maxns)
  menu.compressor=dialog.AddCheckbox("(C) Compressor",ped.run_compressor)
  menu.compressor_loops=dialog.AddTextbox("C loops",ped.compressor_loops)
  menu.quaker=dialog.AddCheckbox("(QR) QuakeR",ped.run_quake)
  menu.quake_loops=dialog.AddTextbox("QR loops",ped.quake_loops)
  menu.voidskiller=dialog.AddCheckbox("(VK) Voids Killer",ped.run_voids_killer)
  menu.voids_killer_loops=dialog.AddTextbox("VK loops",ped.voids_killer_loops)
  menu.lsquake=dialog.AddCheckbox("(LSQ) LS Quake",ped.run_ls_quake)
  menu.ls_quake_loops=dialog.AddTextbox("LSQ loops",ped.ls_quake_loops)
  menu.total_loops=dialog.AddTextbox("Total loops",ped.total_loops)
--menu.relax_option=dialog.AddTextbox("Relax option (0-5)",ped.relax_protein_option)
  menu.relax_option=dialog.AddSlider("Relax option",1,0,5,0)
  menu.order=dialog.AddTextbox("Order",ped.bander_order)
  menu.run=dialog.AddButton("Run",1)
  rc=dialog.Show(menu)
  if rc==1 then
    ped.run_quake_slice=menu.quakeslice.value
    ped.quake_slice_maxns=tonumber(menu.quake_slice_maxns.value)
    ped.run_compressor=menu.compressor.value
    ped.compressor_loops=tonumber(menu.compressor_loops.value)
    ped.run_quake=menu.quaker.value
    ped.quake_loops=tonumber(menu.quake_loops.value)
    ped.run_voids_killer=menu.voidskiller.value
    ped.voids_killer_loops=tonumber(menu.voids_killer_loops.value)
    ped.run_ls_quake=menu.lsquake.value
    ped.ls_quake_loops=tonumber(menu.ls_quake_loops.value)
    ped.total_loops=tonumber(menu.total_loops.value)
    ped.relax_protein_option=tonumber(menu.relax_option.value)
    ped.bander_order=menu.order.value
  end
end

-- ============== functions for quake slice, compressor, quake, voids killer, and LS Quake ============

function FuseEnd()
  behavior.SetClashImportance(.6)
  Wiggle()
  Wiggle("s",1)
  Wiggle()
end

function SaveBest()
  local s=current.GetScore()
  local g=s-bestScore
  if g>0 then
    print(string.format("Gained another %.3f pts.",g))
    bestScore=s
    save.Quicksave(3)
  end
end

function Wiggle(how, iters, minppi) --score conditioned recursive wiggle/shake
  if how==nil then how="wa" end
  if iters==nil then iters=6 end
  if minppi==nil then minppi=0.04 end

  if iters>0 then
    iters=iters-1
    local sp=current.GetScore()
    if how == "s" then structure.ShakeSidechainsSelected(1)
      elseif how == "wb" then structure.WiggleSelected(2,true,false)
      elseif how == "ws" then structure.WiggleSelected(2,false,true)
      elseif how == "wa" then structure.WiggleSelected(2,true,true)
    end
    local ig=current.GetScore()-sp
    if ig > minppi then return Wiggle(how, iters, minppi) end --recurence tail call ;]
  end
end

function bandstr(str) --set all band strengths
  local nbands=band.GetCount()
  local oldbands=ped.bands_at_start
  if nbands>oldbands then
    print('Changing strength of bands '..(1+oldbands)..'-'..nbands..' to '..str..' in bandstr.')
    for i=1+oldbands, nbands do
      band.SetStrength(i,str)
    end -- for i
  end -- if nbands
end

function SaveRB(slot)
  if slot==nil then slot=4 end
  save.Quicksave(slot)
  recentbest.Restore()
  SaveBest()
  save.Quickload(slot)
end


function AllLoop() --turning entire structure to loops
  local seg
  local ok=false
  for i=1, segCnt do
    local s=structure.GetSecondaryStructure(i)
    if s~="L" then
      save.SaveSecondaryStructure()
      ok=true
      break
    end
  end
  if ok then
    selection.SelectAll()
    for seg=1,structure.GetCount() do structure.SetSecondaryStructure(seg,"L") end
  end
end

function qStab()
  behavior.SetClashImportance(0.1)
  Wiggle("s",1)
  selection.SelectAll()
  if fastQstab==false then
    behavior.SetClashImportance(0.4)
    Wiggle("wa",1)
    behavior.SetClashImportance(1)
    Wiggle("s",1)
  end
  behavior.SetClashImportance(1)
  Wiggle()
end

function FuzeEnd()
    behavior.SetClashImportance(1)
    Wiggle("wa",1)
    Wiggle("s",1)
    Wiggle()
    SaveBest()
end

function Fuze1(ci1,ci2)
    behavior.SetClashImportance(ci1)
    Wiggle("s",1)
    behavior.SetClashImportance(ci2)
    Wiggle("wa",1)
end

function Fuze2(ci1,ci2)
  behavior.SetClashImportance(ci1)
  Wiggle("wa",1)
  behavior.SetClashImportance(1)
  Wiggle("wa",1)
  behavior.SetClashImportance(ci2)
  Wiggle("wa",1)
end

function reFuze(scr,slot)
  local s=current.GetScore()
  if s<scr then
    save.Quickload(slot)
  else
    scr=s
    save.Quicksave(slot)
  end
  return scr
end

function Fuze(slot)
  if slot==nil then slot=4 end
  local scr=current.GetScore()
  save.Quicksave(slot)
  selection.SelectAll()
  Fuze1(0.3,0.6) FuzeEnd()
  scr=reFuze(scr,slot)
  Fuze2(0.3,1) SaveBest()
  scr=reFuze(scr,slot)
  Fuze1(0.05,1) SaveBest()
  scr=reFuze(scr,slot)
  Fuze2(0.7,0.5) FuzeEnd()
  scr=reFuze(scr,slot)
  Fuze1(0.07,1) SaveBest()
  reFuze(scr,slot)
end

function MakeBands(num)
  local done=false
  local s1=math.random(segCnt)
  local s2=math.random(segCnt)

  if s1>s2 then s1,s2=s2,s1 end

  local d=structure.GetDistance(s1,s2)
  if s2-s1>=minSkip and d>=minDist then done=true end

  if (done and #useRegions>0) then
    done=false
    for a=1,#useRegions do
      local ss=useRegions[a][1]
      local se=useRegions[a][2]
      if (s1>=ss and s1<=se) or (s2>=ss and s2<=se) then
        done=true
        break
      end
    end
  end

  local ss1=structure.GetSecondaryStructure(s1)
  local ss2=structure.GetSecondaryStructure(s2)

  if (done and noLoops) then
    if ss1=="L" or ss2=="L" then done=false end
  end

  if (done and structure1) then
    if ss1=="L" and ss2=="L" then done=false end
  end

  if (done and S2H) then
    if (ss1=="E" or ss2=="E") and (ss1=="H" or ss2=="H") and (ss1~=ss2) then
      done=true
    else
      done=false
    end
  end

  if done then
    num=num-1
    mkBand1(s1,s2)
    bands[#bands+1]={s1,s2}
  end

  if num>0 then return MakeBands(num) end --tail call
end

function mkBand1(s1,s2)
  band.AddBetweenSegments(s1,s2) 
  -- what if line above fails to add a band?
  local len=structure.GetDistance(s1,s2)
  if decomp==true then
    len=len+compressFrac
  else
    len=len-compressFrac
  end
  if len<3.0 then len=3.0 end
  local nbands=band.GetCount()
  local oldbands=ped.bands_at_start
  if nbands>oldbands then
     band.SetGoalLength(nbands,len)
  else
     print('ERROR in mkBand1('..s1..','..s2..') nbands='..nbands..' <= oldbands='..oldbands..'.')
  end -- if nbands
end

function Repeat_bands(bands)
  MyBandDeleteAll()
  for i=1,#bands do
    local s1=bands[i][1]
    local s2=bands[i][2]
    mkBand1(s1,s2)
  end
end

function Bandstr(str) --set all band strengths
  if str<0.1 then str=0.1 end
  if str>10 then str=10 end
  local nbands=band.GetCount()
  local oldbands=ped.bands_at_start
  if nbands>oldbands then
    print('Changing strength of bands '..(1+oldbands)..'-'..nbands..' to '..str..' in Bandstr.')
    for i=1+oldbands,nbands do
      band.SetStrength(i, str)
    end -- for i
  end -- if nbands
end

function getDist()
  if #distances<1 or distScore~=current.GetScore() then --run only if needed
    print("Calculating distances...")
    for i=1,segCnt do --filling table
      distances[i]={} --need to delclare second dimension
      for j=i+1,segCnt do --not counting from beginning - not need :)
        distances[i][j]=structure.GetDistance(i,j)
      end
    end
    print("Done")
    distScore=current.GetScore()
  end
end

function dist(a,b)
  if a==b then return 0 end
  if a>b then a,b=b,a end
  return distances[a][b]
end

function Sort(tab,items)
  for x=1,items do
    for y=x+1,#tab do
      if tab[x][2]>tab[y][2] then
        tab[x],tab[y]=tab[y],tab[x]
      end
    end
  end
  return tab
end

function mkBand2(a) --make band if found void in area of that segment
  print(string.format("Banding segment %i", a))
  getDist()
  local t={}--there we store possible sehments
  for b=1,segCnt do --test all segments
    local ab=dist(a,b) --distance between segments
    if ab~=nil then -- ab was nil for segment 31 of puzzle 644 11/3/12
    if ab>minLength then --no void if less
      local void=true
      for c=1,segCnt do --searhing that is any segment between them
        local ac=dist(a,c)
        local bc=dist(b,c)
        if ac~=0 and bc~=0 and ac<ab and bc<ab and ac>4 and bc>4 then
          if ac+bc<ab+1.5
            then void=false break --no void there for sure
          end -- if ac
        end -- if ac
      end -- for c
 
      if void==true then
        if math.abs(a-b)>=minDist then
          t[#t+1]={a,b}
        end -- if math
      end -- if void
    end -- if ab>minLength
    end -- if ab~=nil
  end -- for b
  if #t>0 then
    print(string.format("Found %i possible bands across voids",#t))
    for i=1,#t do
      local s=t[i][1]
      local e=t[i][2]
      band.AddBetweenSegments(s,e)
      -- what if line above fails to add a band?
      local d=structure.GetDistance(s,e)
      d=d-compFrac
      if d<3 then d=3 end
      if d>20 then d=20 end
      local nbands=band.GetCount()
      local oldbands=ped.bands_at_start
      if nbands>oldbands then
         band.SetGoalLength(nbands,d)
      else 
         print('ERROR in mkBand2('..a..') nbands='..nbands..' <= oldbands='..oldbands..'.')
      end -- if nbands
    end -- for i
  else
    print("No voids found")
  end -- if #t
end -- mkBand2(a)

function Pull(sn)
  MyBandDeleteAll()
  mkBand2(sn)
  if band.GetCount()>ped.bands_at_start then
    selection.SelectAll()
    behavior.SetClashImportance(pullingCI)
    recentbest.Save()
    local loss=math.abs(math.floor(current.GetScore()*maxLoss/100))
    for str=lastBS,maxBS,0.11 do--search enough band strength to move
      if normal then recentbest.Restore() end --because sometimes it makes points during pull :D
      ss=current.GetScore()
      bandstr(str)
      structure.WiggleSelected(1,true,false)
      if ss-current.GetScore()>loss then
        lastBS=str-0.1
        if lastBS<minBS then lastBS=minBS end
        break
      end
    end
    MyBandDeleteAll()
    SaveRB() --because sometimes it missing fractions
    print("Stabilizing...")
    structure.WiggleSelected(1)
    recentbest.Save() --after pulling
    qStab()
    if bestScore-current.GetScore()<doFuze then
      SaveBest()
      print("Fuzing....")
      Fuze(4)
    end
    SaveBest()
    save.Quickload(3) --load best state
    tmpscore=current.GetScore()
    print(string.format("Current score: %.3f Total gain: %.3f",trunc3(tmpscore),tmpscore-qsc))
    if tmpscore>qsc then listdiff() end
  end
end

function ShellSort ( ids , distances , n )
    -- Adapted from Numerical Recipes in C
    local  inc = 1
    repeat
      inc = inc * 3 + 1
    until inc > n
    repeat
        inc = inc / 3
        inc = inc - inc % 1
        for i = inc + 1 , n do
            v = distances [ i ]
            w = ids [  i ]
            j = i
            flag = false
            while ( flag == false and distances [ j - inc ]  > v ) do
                distances [ j ] = distances [ j - inc ]
                ids [ j ] = ids [ j - inc ]
                j = j - inc
                if ( j <= inc ) then
                  flag = true
                end
            end
          distances [ j ] = v
           ids [ j ] = w
        end
    until inc <= 1
end

function test_for_minimum(i,seg,seg_distances)
  if ((i~=seg) and
     (seg_distances[i] <= seg_distances[i-1]) and
     (seg_distances[i] <= seg_distances[i-2]) and
     (seg_distances[i] <= seg_distances[i+1]) and
     (seg_distances[i] <= seg_distances[i+2])) then
    return true
  else
    return false
  end
end

function test_for_maximum(i,seg,seg_distances)
  if ((seg_distances[i] >= seg_distances[i-1]) and
     (seg_distances[i] >= seg_distances[i-2]) and
     (seg_distances[i] >= seg_distances[i+1]) and
     (seg_distances[i] >= seg_distances[i+2])) then
    return true
  else
    return false
  end
end

function get_band_segment_data(seg)
  seg_distances = {}
  -- Get the function that looks like this:
  -- x = residue number
  -- y = distance of that segment from seg
  -- ignore the degenerate minimum
  n_segments = 0
  segment_ids = {}
  segment_distances = {}
  for i=1,n_residues do
    segment_distances[i] = structure.GetDistance(seg,i)
  end
  for  i = 3,n_residues-3 do
    k=segment_distances[i]
    if ((band_to_minima == true) and
       (test_for_minimum(i,seg,segment_distances) == true )) then
      n_segments=n_segments+1
      segment_ids[n_segments]=i
      segment_distances[n_segments]=segment_distances[i]
    elseif ((band_to_minima == false) and
           (test_for_maximum(i,seg,segment_distances) == true)) then
      n_segments=n_segments+1
      segment_ids[n_segments]=i
      segment_distances[n_segments]=segment_distances[i]
      end
  end
  ShellSort(segment_ids,segment_distances,n_segments)
  if (n_segments >= max_bands_per_residue) then
    n_c=max_bands_per_residue
  else
    n_c=n_segments
  end
  for i=1,n_c do
    if ((segment_distances[i] <= max_segment_length) and
       (seg < min_exclusion or seg > max_exclusion) and
       (segment_ids[i] < min_exclusion or segment_ids[i] > max_exclusion)) then
      n_candidates=n_candidates+1
      candidate_ids_start[n_candidates]=seg
      candidate_ids_end[n_candidates]=segment_ids[i]
      candidate_distances[n_candidates]=segment_distances[i]
    end
  end
end

function create_bands()
  if ( scale_band_strength == true ) then
    total_band_length=0
    for i=1,n_candidates do
      total_band_length=total_band_length + candidate_distances[i]
    end
    average_band_length = total_band_length/n_candidates
  end

  local oldbands=ped.bands_at_start
  for i=1,n_candidates do
    band.AddBetweenSegments(candidate_ids_start[i],candidate_ids_end[i])
    -- what if line above fails to add a band?
    if (scale_band_strength == true) then
      band.SetStrength(i+oldbands,band_strength*average_band_length/candidate_distances[i])
    else
      band.SetStrength(i+oldbands,band_strength)
    end
  end
end

function Quakel()--a Quake
  --do_unfreeze_all()
  save.Quicksave(3)
  print(string.format("Starting Quake by slice. Start score: %.3f",trunc3(qsc)))
  selection.DeselectAll()
  local loss=math.abs(math.floor(current.GetScore()*maxLoss/100))
  print(string.format("Pulling until loss of %i pts.",loss))
  for ns=1,maxns do
    for sl=1 ,ns do
      MyBandDeleteAll()
      maxl=0 unl=0 deuxl=0
      print("Calculating lengths for slice....")
      for i=1,segCount-1 do
        local liun=structure.GetDistance(i,un)
        for j=2,segCount do
          local ljun=structure.GetDistance(j,un)
          local lij=structure.GetDistance(i,j)
          if liun>=(sl-1)*max /ns and liun<=sl*max/ns
             and ljun>=(sl-1)*max /ns and ljun<=sl*max/ns then
            if lij> maxl then
              unl=i deuxl=j maxl=lij
            end -- great than maxl
          end -- in the slice
        end --j
      end --i
      if unl>0 and deuxl>0 and maxl>0 then --check calculations
        for y=1,segCount do
          min=9999999999
          local lyun=structure.GetDistance(y,un)
          if lyun>=(sl-1)*max /ns and lyun<= sl*max/ns and structure.GetDistance(y,unl) + structure.GetDistance(deuxl,y) >maxl*1.1 then
            for x=1,segCount do
              local lxun=structure.GetDistance(x,un)
              if lxun>=(sl-1)*max /ns and lxun<=sl*max/ns and structure.GetDistance(x,unl) + structure.GetDistance(deuxl,x) < maxl*1.1 then
                local lxy=structure.GetDistance(y,x)
                if lxy< min then
                  mins=x min=lxy
                end -- less than min
              end -- x in the slice and in the middle
            end -- x
            if structure.GetDistance(y,mins)>4.5 and math.abs(y-mins)>5 then
              band.AddBetweenSegments(y,mins) --band only if segments are not so close
            end
          end -- y in the slice and in the exterior
        end -- y
      end
      if band.GetCount()>ped.bands_at_start then --check do we have any band
        selection.SelectAll()
        behavior.SetClashImportance(pullingCI)
        recentbest.Save()
        print(string.format("Slice %i/%i/%i min: %i max: %i distance: %.3f",sl,ns,maxns,unl,deuxl,maxl))
        for str=0.1,1.2,0.07 do--search enough band strength to move
          recentbest.Restore()--because sometimes it makes points during pull :D
          ss=current.GetScore()
          bandstr(str)
          structure.WiggleSelected(1,true,false)
          if ss-current.GetScore()>loss then break end
        end
        behavior.SetClashImportance(1)
        MyBandDeleteAll()
        SaveRB() --because sometimes it missing fractions
        recentbest.Save() --after pulling
        qStab()
        if bestScore-current.GetScore()<doFuze then
          SaveBest()
          Fuze()
        end
        SaveBest()
        save.Quickload(3) --load best state
        tmpscore=current.GetScore()
        print(string.format("Current score: %.3f Total gain: %.3f",trunc3(tmpscore),tmpscore-qsc))
        if tmpscore>qsc then listdiff() end
      end
    end -- slice
  end -- number of slice
end -- quake

function Compressor(numTries)
  print(string.format("Starting Compressor v2, %i loops.",numTries))
  if normal==false then print("Using exploration options.") end

  if allLoop then AllLoop() end
  save.Quicksave(3)
  save.Quicksave(7) --save state for losing
  badOnes=0
  recentbest.Save()

  local best=false
  for i=1,numTries do
    loss=current.GetScore()*percLoss/100 --pulling till we lost some points
    MyBandDeleteAll()
    if best==true and repeatGood==true then
      print("Repeating last bands.")
      if modualtor==true then --switch
        decomp= not decomp
      end
      Repeat_bands(bands)
    else
      bands={}
      MakeBands(numBands)
    end -- if best
    local ls=current.GetScore()
    local tx="ompressing."
    if decomp==true then tx="Dec"..tx else tx="C"..tx end
    print(string.format("Loop %i of %i started. %s Current score: %.3f",i,numTries,tx,trunc3(ls)))

    behavior.SetClashImportance(pullingCI)
    selection.SelectAll()
    recentbest.Save()
    if slowBands then
      for str=lastBands,bandStr,0.07 do--search enough band strength to move
        recentbest.Restore() --because sometimes it makes points during pull :D
        local ss=current.GetScore()
        Bandstr(str)
        Wiggle("wb",1)
        if ss-current.GetScore()>loss then
          lastBands=str-0.1
          break
        end
      end -- for str
    else
      Bandstr(bandStr)
      Wiggle("wb",1)
    end -- if slowBands
    SaveRB(4) --because sometimes it missing fractions
    MyBandDeleteAll()
    qStab()
    if useFuze and current.GetScore()> (bestScore-fuzeScore) then
      SaveBest()
      Fuze(4)
    else
      SaveBest()
    end -- if useFuze

    if maxLoss>0 then
      if current.GetScore()>bestScore-maxLoss then
        save.Quicksave(7) --save actual "best" solution
      else
        badOnes=badOnes+1
        if badOnes >=reloadBest then
          badOnes=0
          save.Quickload(3)
        else
          save.Quickload(7) --load "second best"
        end -- if badOnes
      end -- if current.GetScore
    else
      save.Quickload(3)
    end -- if maxLoss
    local es=current.GetScore()
    if es>ls then best=true else best=false end --repeating all bands making points
    print(string.format("Loop %i gain %.3f ;total gain %.3f",i,es-ls,bestScore-sscore))
    if es>ls then listdiff() end
    if modualtor==true then --switch
      decomp= not decomp
    end
  end -- for i
  save.Quickload(3)
  if allLoop then save.LoadSecondaryStructure() end
  print(string.format("Total gain: %.3f",current.GetScore()-sscore))
end

function QuakeR()
  lastBS=minBS
  --do_unfreeze_all()
  local qsc=current.GetScore()
  save.Quicksave(3)
  print(string.format("Starting QuakeR v3 %i passes. Start score: %.3f",qLoops,trunc3(qsc)))
  for xx=1,qLoops do
    print(string.format("Pass %i of %i",xx,qLoops))
    local start=math.floor(math.random(segCnt/10)) --first band somewhere in first 10% of protein
    local len=math.floor(math.random(segCnt/2-5))+10
    local step=math.floor(math.random(segCnt/2-5))+10
    local loss=math.abs(math.floor(current.GetScore()*maxLoss/100))
    MyBandDeleteAll()
    print(string.format("Bands from sgmnt %i of len %i every %i sgmnt.",start,len,step))
    print(string.format("Pulling until loss of %i pts.",loss))
    for x=start,segCnt,step do
      for y=start+len,segCnt,step do
        if y<=segCnt then band.AddBetweenSegments(x,y) end
      end -- for y
    end -- for x
    selection.SelectAll()
    behavior.SetClashImportance(pullingCI)
    recentbest.Save()
    for str=lastBS,maxBS,0.07 do--search enough band strength to move
      if normal then recentbest.Restore() end--because sometimes it makes points during pull :D
      ss=current.GetScore()
      bandstr(str)
      structure.WiggleSelected(1,true,false)
      if ss-current.GetScore()>loss then
        lastBS=str-0.1
        if lastBS<minBS then lastBS=minBS end
        break
      end -- if ss
    end -- for str
    MyBandDeleteAll()
    SaveRB() --because sometimes it missing fractions
    print("Stabilizing...")
    behavior.SetClashImportance(1)
    recentbest.Save() --after pulling
    qStab()
    if bestScore-current.GetScore()<doFuze then
      SaveBest()
      print("Fuzing....")
      Fuze(4)
    end -- if bestScore
    SaveBest()
    save.Quickload(3) --load best state
    tmpscore=current.GetScore()
    print(string.format("Current score: %.3f Total gain: %.3f",trunc3(tmpscore),tmpscore-qsc))
    if tmpscore>qsc then listdiff() end
  end -- for xx
  print(string.format("Total QuakeR gain: %.3f",current.GetScore()-qsc))
end

function KillVoids()
  if se==nil then se=segCnt end
  save.Quicksave(3)
  qsc=current.GetScore()
  lastBS=minBS
  local s={}
  print(string.format("Starting Void Killer v0.4 by rav3n_pl Start score: %.3f",trunc3(qsc)))
  getDist() --calculate distances now
  se=math.min(se,structure.GetCount())

  for i=ss,se do
    s[#s+1]={i,current.GetSegmentEnergyScore(i)}
  end
  if Vworst==true then
    Sort(s,#s)
  end
  if Vrandom==true then
    for i=1,#s do
      local r=math.random(#s)
      s[i],s[r]=s[r],s[i]
    end
  end

  for i=1,#s do
    print(string.format("Kill %i of %i started.",i,#s))
    Pull(s[i][1])
  end
end

function Quake(seg)
  selection.SelectAll()
  init_score=current.GetScore()
  structure.WiggleSelected(1,true,false)
  score_after_wiggle =current.GetScore()
  score_reduction_after_wiggle=init_score - score_after_wiggle
  if (monit) then
    print(string.format("seg %i score loss %.3f",seg,score_reduction_after_wiggle))
  end
  
  -- Adaptively change band strength based on how much the previous wiggle overshot/undershot the target. Could
  -- be made much more elaborate.
  
  if (score_reduction_after_wiggle>target_wiggle_reduction*1.2) then
    band_strength=band_strength*(0.85+0.1*math.random())
  elseif (score_reduction_after_wiggle<target_wiggle_reduction*0.8 ) then
    band_strength=band_strength+(1-band_strength)*(0.05+0.1*math.random())
  end
  if (monit) then
    print(string.format("band_strength %.3f",band_strength))
  end
  MyBandDeleteAll()
  structure.ShakeSidechainsSelected(1)
  gas=current.GetScore()-score_after_wiggle
  if (monit) then
    print(string.format("      gas %.3f", gas ))
  end
  structure.WiggleSelected(8)
  if (do_second_shake==true) then
    score_a2w=current.GetScore()
    structure.ShakeSidechainsSelected(1)
    score_a2s=current.GetScore()
    if (score_a2s-score_a2w>1) then
      structure.WiggleSelected(4)
    end
 end
end
-- ============== End of functions ===============

ped.note_number=structure.GetCount()
for seg=structure.GetCount(),1,-1 do
  if structure.GetNote(seg)~="" then break end
  ped.note_number=seg
end
print(string.format("Recording MB results in Note for segment %i",ped.note_number))
ped.starting_score=current.GetScore()
structure.SetNote(ped.note_number,string.format("(%s) %.3f + MB",user.GetPlayerName(),trunc3(ped.starting_score)))

if ped.use_gui then multibander_gui() end

xref={}
xref["QS"]="Quake Slice"
xref["C"]="Compressor"
xref["QR"]="QuakeR"
xref["VK"]="Voids Killer"
xref["LSQ"]="LS Quake"
ped.bander_order1={}
list=ped.bander_order..","
for w in string.gmatch(list,"%a+") do
  if xref[w]~=nil then
    ped.bander_order1[#ped.bander_order1+1]=xref[w]
  end
end

segCount=structure.GetCount() --always the same

bestScore=current.GetScore()

start_clock=os.clock()
start_time=os.time()
report_options()
record_score("Initial")
ped.disable={}
ped.disable["Quake Slice"]=false
ped.disable["Compressor"]=false
ped.disable["QuakeR"]=false
ped.disable["Voids Killer"]=false
ped.disable["LS Quake"]=false
ped.last_score={}
ped.best_score=-99999999
ped.protein_was_relaxed=false
save.Quicksave(50)

for pedloop=1,ped.total_loops do

  ped.loop_count=pedloop
  ped.loop_start_score=current.GetScore()

  for bander_count=1,#ped.bander_order1 do
    ped.bander=ped.bander_order1[bander_count]
    if ped.run_quake_slice and ped.bander=="Quake Slice" then make_a_decision("Quake Slice") end
    if ped.run_quake_slice and ped.bander=="Quake Slice" and not ped.disable["Quake Slice"] then

      print("============= Quake Slice v1.3 =============")

      ped.score1=current.GetScore()
      max=0 un=0 deux=0
      print("Searching for central axe...")
      for i=1,segCount-1 do
          for j=2,segCount do
             local len=structure.GetDistance(i,j)
             if len > max then
                 un=i
                 deux=j
                 max=len
             end
          end
      end

      maxLoss=0.5 --minimum percentage loss when pulling. ie 5% of 10000 is 500pts
      pullingCI=0.9 --clash importance while pulling
      doFuze=10 --run fuze when score after qstabilize is close to best
      maxns=ped.quake_slice_maxns -- maximum number of slices
      fastQstab=false --if true ony 1s1w as stabilize after pull


      if ped.quake_slice_un then
        qsc=current.GetScore() --starting score -- moved here for code level 3.03 V2
        bestScore=current.GetScore() -- added for multibander v2.04 -- move here for code level 3.03 V2
        Quakel()
      end

      if ped.quake_slice_deux then
        qsc=current.GetScore() --starting score -- moved here for code level 3.03 V2
        bestScore=current.GetScore() -- added for multibander v2.04 -- move here for code level 3.03 V2
        un=deux --from the other side now
        Quakel()
      end

      record_score("Quake Slice")
      ped.gains["Quake Slice"]=current.GetScore()-ped.score1
      ped.last_score["Quake Slice"]=current.GetScore()

    end


    if ped.run_compressor and ped.bander=="Compressor" then make_a_decision("Compressor") end
    if ped.run_compressor and ped.bander=="Compressor"  and not ped.disable["Compressor"] then

      print("============= Rav3n_pl Compressor v2.2 =============")

      ped.score1=current.GetScore()
      normal=true --set false for exploration puzles
--rest of options at end

      segCnt=structure.GetCount()
      while structure.GetSecondaryStructure(segCnt)=="M" do segCnt=segCnt-1 end
      bands={}

      sscore=current.GetScore() --starting score

----------- options below VVVVV ---------------------

      useRegions={  --set regions that have to be used in all bands
      --{1,20}
      --{25,55}
      --{100,210}
      }

      compressFrac=5 --make bands shorter by that much points
      minSkip=15 --minimum segment distance between banded segments
      numBands=8 --how many bands use at once

      fastQstab=true --true --only 1 shake and 1 wiglle when true

      useFuze=true -- use Fuze after qStab
      fuzeScore=-1 --how close we have to be to run Fuze

      allLoop=ped.compressor_allLoop --work in all-loop mode. sometimes work better than structure mode :)

      S2H=false --all bands are between sheets and helixes
      structure1=false --all bands have at least one end on structure (not loop) -- modified name from structure to structure1 for foldit V2
      noLoops=false --band can`t have any end on loop

      lastBands=0.3 --starting band str
      bandStr=1.0 -- max band strength
      minDist=7 --mimum band length
      slowBands=true --bands strength are raised few times till score drops by %

      percLoss=1 --pulling stops when score drops by that percent (ie 2=200pts when 10k pts)
      repeatGood=false --repeating good bands

      pullingCI=1.0 --clash importance during pull
      maxLoss=30 --maximum acceptable LOSS after  pull/stabilize
      reloadBest=5 --reload best solution after that many worst ones

      decomp=false --true --true for DeCompression instead of compression
      modualtor=ped.compressor_modulator --true --changing compression/decompression each loop

      pd_best=-99999
      save.Quicksave(6) -- store starting solution
      slots={}
      slots={3,7} -- slots 3,7 used by compressor

      bestScore=current.GetScore() -- (need to add this in loop to get compressor to accumulate the best scores correctly)
      sscore=current.GetScore() --starting score (need this in the loop for compressor to score correctly)
      Compressor(ped.compressor_loops)

      record_score("Compressor("..ped.compressor_loops..")")
      ped.gains["Compressor"]=current.GetScore()-ped.score1
      ped.last_score["Compressor"]=current.GetScore()

    end

    if ped.run_quake and ped.bander=="QuakeR" then make_a_decision("QuakeR") end
    if ped.run_quake and ped.bander=="QuakeR" and not ped.disable["QuakeR"] then

      print("============= QuakeR =============")


      ped.score1=current.GetScore()
      normal=true --set false for exploration puzzles
      -- rest of Options at end


      segCnt=structure.GetCount()
      while structure.GetSecondaryStructure(segCnt)=="M" do segCnt=segCnt-1 end

      --- VVVVVVVVVV OPTIONS

      minBS=0.3 --starning minimum bands strength
      maxBS=1.0 --maximum bands strength
      maxLoss=1 --minimum percentage loss when pulling. ie 5% of 10000 is 500pts
      pullingCI=1.0 --clash importance while pulling
      fastQstab=false --true for 1s1w as stabilize (faster)
      doFuze=10 --run fuze when score after qstabilize is close to best. If negative  run fuze only if  gain after qstab

      bestScore=current.GetScore() -- (need to add this in loop to get compressor to accumulate the best scores correctly)
      qLoops=ped.quake_loops  -- <<<OVER THERE! do more! much more!
      QuakeR()

      record_score("QuakeR("..ped.quake_loops..")")
      ped.gains["QuakeR"]=current.GetScore()-ped.score1
      ped.last_score["QuakeR"]=current.GetScore()

    end

-- =======================================

    if ped.run_voids_killer and ped.bander=="Voids Killer" then make_a_decision("Voids Killer") end
    if ped.run_voids_killer and ped.bander=="Voids Killer" and not ped.disable["Voids Killer"] then

      print("============= Voids Killer v0.4 =============")


      ped.score1=current.GetScore()
      normal=true --set false if want use exploration score not energy score
      --rest of options at end

      segCnt=structure.GetCount()
      while structure.GetSecondaryStructure(segCnt)=="M" do segCnt=segCnt-1 end

      bestScore=current.GetScore()
      distances={}
      distScore=current.GetScore()


      --- OPTIONS THERE! vvv

      minDist=15 --minimum segment distance
      minLength=10 --minimum spatial distance

      Vworst=true --start form worst scoring one
      Vrandom=false --true --go it in random order
      -- if both are false then it will band form 1st to last

      doFuze=-1 --run fuzes when close to saved best (negative - run if fain after qstab)
      pullingCI=0.8 --clash improtance during pull
      fastQstab = true --1s1w then true, double if false

      minBS=0.3 --startng band strength
      maxBS=1 --maximum bands strength
      compFrac=3 --how much shorter band should be
      maxLoss=1 --pull till perc loss (ie 1 on 10`000 ptsy = pull till 100 pt loss after pull)

      ss=1 --from segment 1
      se=nil --to last segment if nil
      se=ped.voids_killer_loops -- set number of loops

      --END OF OPTOPNS ^^^

      KillVoids()

      record_score("Voids Killer")
      ped.gains["Voids Killer"]=current.GetScore()-ped.score1
      ped.last_score["Voids Killer"]=current.GetScore()

    end

-- =======================================
    if ped.run_ls_quake and ped.bander=="LS Quake" then make_a_decision("LS Quake") end
    if ped.run_ls_quake and ped.bander=="LS Quake" and not ped.disable["LS Quake"] then

      print("============= LS Quake =============")
      ped.score1=current.GetScore()

      use_ranked_score = 1                       -- if 1, use the score being used for ranking (best normally). If 0, use the
                                                                     --  old style scoring system
      max_bands_per_residue = 6           -- Increase this for more bands at each segment
      max_segment_length = 999               -- Don't create bands greater than this length. Set to 999 or something if you're
                                                                     -- not worried about the potential for long bands creeping in,
      min_exclusion =  999                       -- These can be set to establish an "exclusion zone", within which no bands will
                                                                   --  be created. This might be used to prevent banding to a wayward loop.
      max_exclusion =  1
      target_wiggle_reduction = 150        -- On wiggle, try and reduce the score by this amount after creating the bands.
                                                                      -- If it's wildly out, reduce/increase band strength
                                                                      -- for the next try (don't try to adjust for the current scramble).
      max_segment_length = 10               -- Don't create bands greater than this length. Set to 999 or something if you're
                                                                     -- not worried about the potential for long bands creeping in,
      n_quakes  = ped.ls_quake_loops
      band_to_minima = true                -- Band to minima or maxima in the distance. See bands_from_segment for meaning.
      do_second_shake = false           -- After creating the bands and pulling with them do a Shake/Wiggle/Shake
                                                                -- if do_second_shake is false (best I think) don't do the last shake.
                                                                -- Speeds things up quite a bit
      scale_band_strength = true        -- If true, band strength will be scaled inversely with length
      monit = false


---------- A few more globals

      candidate_ids_start = {}
      candidate_ids_end = {}
      candidate_distances = {}
      n_candidates = 0
      best_score = 0

      n_residues = structure.GetCount()
      band_strength = 0.5
      behavior.SetClashImportance(1)

      recentbest.Save()

      best_score = current.GetScore()

      print(string.format(" Initial score: %.3f",trunc3(best_score)))

      idx = 1
      gel=1

      for i =  1, n_quakes  do

        recentbest.Restore()

        -- Need to record any improvements from the previous iteration
        score = current.GetScore()
        if ( score > best_score ) then
           gain=score-best_score
           best_score = score
           print(string.format("    Gain: %.3f Score: %.3f",gain,trunc3(score)))
           listdiff()
        end

        freeze.UnfreezeAll()
        MyBandDeleteAll()
        n_candidates = 0
        debut=0
        if (idx==1) then
          debut=1
        else
          if ( structure.GetSecondaryStructure(idx)~=structure.GetSecondaryStructure(idx-1) ) then
            debut=idx
          end
        end -- if idx

        if ( debut>0) then
          while (idx < n_residues+1 and structure.GetSecondaryStructure(idx)==structure.GetSecondaryStructure(debut)) do
            get_band_segment_data  ( idx )
            if gel ==1 then
              selection.DeselectAll()
              selection.Select(idx)
              freeze.FreezeSelected(true,false)
            end -- if gel

            idx=idx+1
          end -- while idx
        end -- if debut
        selection.SelectAll()
        if ( n_candidates > 0 ) then
          create_bands ()
          sc=current.GetScore()
          print(string.format("  Iteration: %i/%i Segments: %i-%i Score: %.3f",i,n_quakes,debut,idx-1,trunc3(sc)))
          Quake ( idx )
        end -- if n_candidates

        if idx>n_residues then
          idx=1
          if gel==1 then gel=0 else gel=1 end
        end -- if idx
      end -- for i
      recentbest.Restore()
      selection.SelectAll()
      freeze.UnfreezeAll()

      record_score("LS Quake")
      ped.gains["LS Quake"]=current.GetScore()-ped.score1
      ped.last_score["LS Quake"]=current.GetScore()

    end

-- =======================================
    if ped.enable_mutate then
      score1=current.GetScore()
      oldseq=getseq() -- get aa sequence that goes with score1
      repeat
        score2=current.GetScore()
        structure.MutateSidechainsAll(1)
        score3=current.GetScore()
      until score3<score2+0.01
      score4=current.GetScore()
      gain=score4-score1
      print(string.format("Mutated from %.3f to %.3f (gain=%.3f)",trunc3(score1),trunc3(score4),gain))
      if gain==gain then -- if gain>0 then -- can successful mutations give gain<=0 ?
         newseq=getseq() -- get aa sequence that goes with score4
         listseqdiffs(oldseq,newseq) -- list the differences between old and new aa sequences
      end -- if gain
    end
-- =======================================
    score=current.GetScore()
    if score>ped.best_score then
      save.Quicksave(50)
      ped.best_score=score
    end
-- =======================================
    if bander_count~=#ped.bander_order1 then
      structure.SetNote(ped.note_number,string.format("(%s) %.3f + MB(%i+) %.3f",user.GetPlayerName(),trunc3(ped.starting_score),pedloop-1,trunc3(current.GetScore())))
    else
      structure.SetNote(ped.note_number,string.format("(%s) %.3f + MB(%i) %.3f",user.GetPlayerName(),trunc3(ped.starting_score),pedloop,trunc3(current.GetScore())))
    end
    save.Quicksave(1)
    save.Quickload(50)
    if bander_count~=#ped.bander_order1 then
      structure.SetNote(ped.note_number,string.format("(%s) %.3f + MB(%i+) %.3f",user.GetPlayerName(),trunc3(ped.starting_score),pedloop-1,trunc3(current.GetScore())))
    else
      structure.SetNote(ped.note_number,string.format("(%s) %.3f + MB(%i) %.3f",user.GetPlayerName(),trunc3(ped.starting_score),pedloop,trunc3(current.GetScore())))
    end
    save.Quicksave(50)
    save.Quickload(1)
-- =======================================


    if bander_count==#ped.bander_order1 then -- last bander has been run
      ped.loop_end_score=current.GetScore()
      loop_gain=ped.loop_end_score-ped.loop_start_score
      print(string.format("Loop gain: %.3f",loop_gain))
      if ped.loop_end_score>ped.best_score then
        save.Quicksave(50)
        ped.best_score=ped.loop_end_score
      end

      ped.protein_was_relaxed=false
      if ped.relax_protein_option > 0 and loop_gain < ped.minimum_loop_gain then
        if ped.relax_from_best then save.Quickload(50) end
        if ped.relax_protein_option == 1 then
          behavior.SetClashImportance(ped.bcis[ped.bci_index])
          structure.WiggleSelected(5,true,false)
          behavior.SetClashImportance(1.0)
          record_score("Relax protein(wiggle down "..ped.bcis[ped.bci_index]..")")
          ped.bci_index=ped.bci_index+1
          if ped.bci_index > #ped.bcis then ped.bci_index=1 end
        elseif ped.relax_protein_option == 2 then
          behavior.SetClashImportance(ped.bcis[ped.bci_index])
          structure.ShakeSidechainsSelected(1)
          behavior.SetClashImportance(1.0)
          record_score("Relax protein(shake down "..ped.bcis[ped.bci_index]..")")
          ped.bci_index=ped.bci_index+1
          if ped.bci_index > #ped.bcis then ped.bci_index=1 end
        elseif ped.relax_protein_option==3 then
          array1=get_segment_scores_array()
          sort_column=2
          table.sort(array1,sort_by_column)
          init_score=current.GetScore()
          for arix=1,#array1 do
            seg=array1[arix][1]
            print(string.format("...rebuilding segment %i with segment score of %.3f",seg,array1[arix][2]))
            selection.DeselectAll()
            selection.Select(seg)
            score=current.GetScore()
            n=1
            repeat
              structure.RebuildSelected(n)
              n=n+1
            until (score~=current.GetScore() or n>10)
            if init_score-current.GetScore() > ped.relaxation_minimum_amount then break end
          end
          selection.SelectAll()
          record_score("Relax protein (rebuild worst)")
        elseif ped.relax_protein_option==4 then
          push_and_pull()
          selection.SelectAll()
          record_score("Relax protein (push and pull)")
        elseif ped.relax_protein_option==5 then
          selection.SelectAll()
          score1=current.GetScore()
          n=1
          repeat
            structure.RebuildSelected(n)
            score2=current.GetScore()
            n=n+1
          until score2~=score1
          repeat
            score1=current.GetScore()
            structure.ShakeSidechainsAll(1)
            score2=current.GetScore()
          until score2<score1+1
          repeat
            score1=current.GetScore()
            structure.WiggleAll(2,false,true)
            score2=current.GetScore()
          until score2<score1+1
          repeat
            score1=current.GetScore()
            structure.WiggleAll(2,true,true)
            score2=current.GetScore()
          until score2<score1+1
          record_score("Relax protein (rebuild whole protein)")
        end
        ped.protein_was_relaxed=true
      end
    end
  end
end

save.Quickload(50)