Recipe: Simulated Annealing 4.1

created by KarenCH

Profile

• Name: Simulated Annealing 4.1
• ID: 103607
• Shared with: Public
• Parent: None
• Children: None
• Created on: June 15, 2020 at 22:19 PM UTC
• Updated on: June 15, 2020 at 22:19 PM UTC
• Description:

  Lots of actions controlled by a Simulated Annealing algorithm for choosing which are to be accepted

Code

--------------------------------------------------------------------
--           Code for Simulated Annealing algorithm
--------------------------------------------------------------------
-- Performs "modifications" of a protein and accepts or rejects them depending 
-- on "temperature" and on how much they improved the score (if + improvement, yes; 
-- if not too bad score drop and temperature is "high", then maybe yes). 
-- Temperature starts high, but drops gradually to 0 by end of run.
--
-- TODO: 
-- 1. Cysteine: 
--     Include band-two-cysteines together (locate closest ones, match together, etc). 
--     Consider throwing out actions that break connections.    
-- 2. Rebuild/remix: 
--     Rebuild with compressor bands.
--     Offer remix
-- 3. Idealize: 
--     Compare with cuts versus without, see if lengths should vary. 
--     Maybe put compressor bands instead of zlb?
-- 4. Contactmaps: 
--     Try more actions such as "put band(s)" + rebuild at low ci; use ZLB-wiggler. 
--     Try a "Til contact is made" banded wiggler. 
--     TEST comp-band mods: are they working propery?
-- 5. LowCI:
--     Should it include extremely lowCI stuff? (mainly idealizing with ZLBs at 0.05)
-- 6. Banders:
--     Experiment more with expanding the protein before actions. rebuild? late stage?
-- 7. Dialog boxes:
--     Missing mutate-related options - particularly a no-mutate list, maybe PROB_MUTATESIMILAR
--     Missing individual actions: wiggles, banders, rebuilders, idealizers (really want this???)
--     USE dlg.foo = dialog.AddTextbox( string label, string value ) 
--         for no-mutate list from user or cysteine-banders
--         Parse value to learn what user wants.
--
-- Code ideas and/or thefts come from: 
--   (tlaloc) tlaloc Random Tug 4.00
--   (spvincent) Helix Twister 1.0, 
--               Loop rebuild 5.0 (NC), 
--               Local Quake 1.0, 
--               Quaking Rebuild V2 1.0
--   (susume2) Cut and Wiggle Everything v0.1
--   (MurloW) Shock v0.1, Fracture v1.6
--   (MurloW, drjr) Idealize by 3+
--   (rav3n_pl) Rav3n_pl GAB BiS v2.0.1 w30 [NC], 
--              Rav3n_pl Voids Killer v0.5 NC, 
--              Rav3n_pl Push v3, LS Quake v1.1
--   (tlaloc, rav3n_pl, Seagat2011, thom001) ST - Glycine Hinge REPOST
--   (drjr) drjr - RandomizeR 1.20, Idealize by 3
--   (BitSpawn) general "how to do low CI" (any failures to do right are mine, not his)
--
--------------------------------------------------------------------------------


----------------------------------------------------------------------
--
---- SIMULATED-ANNEALING-SPECIFIC STUFF
--
----------------------------------------------------------------------

-- Simulated Annealing "core" values
-- The three main phases each perform 2 SA runs.
SA_NUMSTARTS_EARLY     = 100     -- no expander
SA_STEPSPERRUN_EARLY   = 15 
SA_NUMSTARTS_EARLY2    = 70      -- some expander
SA_STEPSPERRUN_EARLY2  = 15 

SA_NUMSTARTS_REBUILD   = 100
SA_STEPSPERRUN_REBUILD = 15 
SA_NUMSTARTS_REBUILD2  = 75 
SA_STEPSPERRUN_REBUILD2= 20     -- try for longer session 2nd go-round

SA_NUMSTARTS_MID       = 1      -- long battle for a few points.
SA_STEPSPERRUN_MID     = 300 
SA_NUMSTARTS_MID2      = 1      -- ugly struggle for a few points.
SA_STEPSPERRUN_MID2    = 300

-- when and how often to pre-expand protein before acting on it
ALLOW_EXPAND_FOR_ACTION   = true
PROB_USE_EXPANDER_EARLY2  = 0.25 -- find out what is best
PROB_USE_EXPANDER_REBUILD = 0.25 -- find out what is best
PROB_USE_EXPANDER_MID1    = 0.5  -- find out what is best
PROB_USE_EXPANDER_MID2    = 1.0  -- seems always good here
PROB_USE_EXPANDER_LATE    = 1.0  -- seems always good here
PROB_USE_EXPANDER = PROB_USE_EXPANDER_EARLY2  -- will get set to action's particular choice

-- SA-Midgame, how often should it target "bad" indices?
PROB_TARGETED = 0.25

-- values for mutate 
PROB_MUTATESIMILAR = 0.25 -- was 0.75


------------ LESS INTERESTING SA STUFF ------------
-- constants for weird SA operations 
SA_STEPSPERRUN_LOW20   = 50     -- CI <= 0.20 -- first half of "LowCI"
SA_NUMSTARTS_LOW20     = 4
SA_STEPSPERRUN_LOW50   = 50     -- CI <= 0.50 -- second half of "LowCI"
SA_NUMSTARTS_LOW50     = 4
SA_STEPSPERRUN_CMAP    = 50     -- special Contact map operations
SA_NUMSTARTS_CMAP      = 4
SA_STEPSPERRUN_LATE    = 50     -- Does fuse. rarely useful
SA_NUMSTARTS_LATE      = 1  

-- values for controlling behavior of fusing during Late stage
ALLOW_FUSE = true                 -- dialog boxes set false for all but OperationChosen=="SA-Late"
FUSEWIGGLE_ITERS   = 2
FUSEWIGGLE_GAIN    = 0.10
PROB_FUSEWIGGLE_CUTFUSE = 0.70
PROB_FUSEWIGGLE_LOCAL = 0.00      -- very slow, very poor success rate
PROB_FUSEWIGGLE_BLUEFUSE = 0.30

-- when to bail if a protein is mostly locked and therefore rebuilds always fail
MAX_REBUILD_FAILS = 20
rebuild_fails = 0

PERFORM_CYSTEINE_BANDING = false
CYSTEINE_BAND_STRENGTH = 5.0    -- 9.0 seems dicy; 5.0 has worked just fine in past
MAX_RI_CYSTEINE_SEPARATION = 15.0  -- if more, then pose is nuts. fail the action that did this
COUNT_CYSTEINE_BANDS_USED    = 3    -- pick 1..4 or break your script

-- if a longish run is being done with non-Zero temperature function, this is a way to bail out
SA_MINGAIN             = 5.0 -- the amount of gain that must be exceeded to not bail out
SA_FAILIFNOGAIN        = 25  -- the step to quit at

--------------------------------------------------------
-- 
----------------    BOILERPLATE    ---------------------
-- 
--------------------------------------------------------

------------------------------------------------------------------------------------------------
-- Standard variables/constants used in my foldit scripts. Some are intended to be changed 
-- as part of setup of scripts, others are user settable, yet more are unlikely to be changed.
-- A separate group, indicated as such, are genuinely intended as constants or as program-
-- controlled variables.
------------------------------------------------------------------------------------------------

-- some oddballs that users might want to set
USENORMALSCORE    = true  -- exploration puzzles would set false
DEBUGRUN          = true  -- mostly goes with DebugPrint, but could get more use later
SHOW_STATS_REPORT = true

-- what to do with user state (bands, cuts) when running, and what cleanup is best at end
PRESERVE_USER_BANDS = true
PRESERVE_USER_CUTS = false
RESTORE_BEST_AT_END = true
PERFORM_CLEANUP_AT_END = true
KILLING_FILTERS = false
EXPAND_BEFORE_ACTION = false

-- simple bander options
BAND_STRENGTH_DEFAULT = 1.0
BAND_STRENGTH_REGION = 1.0
BAND_MINCOMPRESS = 0.80
BAND_MINEXPAND = 1.05
BAND_MAXEXPAND = 1.20
BANDEDWIGGLE_TARGET = 500.0
BIS_LENGTH = 5.0                    -- max length of a normal BIS
BIS_LENGTH_LONG = 8.0               -- max length of a "long" BIS
BAND_TO_SIDECHAINS = true
BAND_ADD_RANDOMS = true
RANDOMBAND_COUNT = 3  -- when random bands are added, this is a plausible count

-- default time to use in shaking
SHAKESIDECHAIN_ITERS   = 1

-- values for controlling in-step wiggle timing and gain
USEITERATIVELONGWIGGLE = FALSE
LONGWIGGLE_ITERS = 12       -- was 10. 14 too slow. consider 10 for early, 12 for rebuild, 14 for midgame?
QUICKWIGGLE_ITERS = 2
QUICKWIGGLE_GAIN  = 0.50    -- for iterative wiggles that happen in quickwiggle time
PROB_QSTAB_EXTRA = 0.50
PROB_LOCALWIGGLE = 1.00 -- I currently believe "before" wiggles should always prefer LocalWiggle over Wiggle.
USE_LOWCI_WIGGLES = false  -- set true by default in low-ci scripts

-- values for controlling how mutations are handled (under-tested)
ALLOW_SPOT_MUTATION = false        -- this is for performing local "spot" changes
MUTATE_NOT_SHAKE = false  -- this is for choosing between shake and mutate

-- values associated with idealize
IDEALIZE_USE_CUTS = true
IDEALIZE_FORBID_UNHEALED_CUTS = true
IDEALIZE_MAXTRIES = 5  -- a "try" fails if cuts not healed. increase ZLBstrength, try again.
IDEALIZE_ZLB_STRENGTH = 2.0
IDEALIZE_TUBE_RADIUS = 4.0
IDEALIZE_TUBE_BIG_RADIUS = 10.0
IDEALIZE_MAX_SEGLEN = 5

-- values associated with rebuilding
REBUILD_FAILSBEFOREQUIT = 75
REBUILD_TRYCOUNT = 5
REBUILD_ADDIDEALIZE = true
REBUILD_USECUTS = true
REBUILD_MAX_SEGLEN = 10          -- prefer even!
REBUILD_MAX_EARLY_SCOREDROP = 3000.0
REBUILD_MAX_MEDIUM_SCOREDROP = 500.0
REBUILD_MAX_LATE_SCOREDROP = 100.0
REBUILD_INNER_TUBE_RADIUS  = 12.0  -- 9.0 has been seen, 7 or 8 recommended for contact map
REBUILD_OUTER_TUBE_RADIUS = 18.0 
REBUILD_CM_INNER_TUBE_RADIUS  = 8.0
REBUILD_CM_OUTER_TUBE_RADIUS = 16.0 

-- Clash Importance
SCALE_MAXCI = 1.0
COREPULLING_CI = 0.70
REGIONMOD_CI = 0.10

-- values for contact map reading 
CONTACTMAP_THRESHOLD = 0.25
CONTACTMAP_TAKEDROPS = false     -- allow drops if they improve contact map
CONTACTMAP_MAX_DROP = 25.0
CONTACTMAP_CONTACTING_GOAL = 0.90
CONTACTMAP_NONCONTACTING_GOAL = 0.70
CONTACTMAP_CONTACTING_STRENGTH = 0.50
CONTACTMAP_NONCONTACTING_STRENGTH = 1.50
MAX_ALLOWED_NONCONTACT = 10.0   -- longer means disbelieve it's a valid contact

-- TailWhip constant
MAX_IDX_FROM_END = 5

----------------------------------------------------------------------------------
------------ VALUES USERS SHOULD NOT BE TOUCHING ---------------------------------
----------------------------------------------------------------------------------

-------------- CONSTANTS THAT PROBABLY SHOULD NEVER CHANGE --------------
-- atoms in an amino acid
BETA_CARBON = 5
TERMINAL_BETA = 6
CENTER_CARBON = 2 --this is the default atom for bands to attach to

-- ideal distances
HELIX_DISTANCE_SKIP4         = 6.3  -- for skips of 4
HELIX_DISTANCE_SKIP5         = 8.6  -- for skips of 5
DISTANCE_SHEET_TO_SHEET      = 4.8
DISTANCE_HELIX_TO_THING      = 8.2
DISTANCE_SHEET_SEG_TO_SEG    = 6.8
DISTANCE_INSHEET_N_TO_NPLUS2 = 9.5
DISTANCE_CYSTEINE_BOND       = 2.05
SCALE_CYSTEINE_OFF_1		 = 1.8

-- QUICKSAVE slots used herein
QS_Start      = 1
QS_CmBest     = 2
QS_Best       = 3
QS_Alternate1 = 4     -- these are for scripts that offer alternate poses
QS_Alternate2 = 5
QS_Alternate3 = 6
QS_Alternate4 = 7
QS_Alternate5 = 8

-- a quickslot stack for any system that wants such (some rebuilder/remixer scripts want)
QS_Base          = 9       -- the base of a stack of slots to operate with
QS_Top           = QS_Base -- the current top of the stack
QS_MAX           = 49      -- the max top point for the stack

QS_RemixSlot1         = 50    -- surely nobody will remix more than 20 slots!
QS_RemixSlotLast      = 69
QS_MultiStartSave1    = 70    -- surely nobody will multistart more than 20 times!
QS_MultiStartSaveLast = 89
QS_Swapslot      = 91 -- if sorting of quickslots is needed, use for swapping
QS_Swapslot2     = 92 -- if sorting of quickslots is needed, use for swapping
QS_ShortUseTemp1 = 93 -- used by low-level actions that stash state for short time
QS_ShortUseTemp2 = 94 -- used by low-level actions that stash state for short time
-- following slots are not safe to borrow if MultiTry or SA are in play
QS_NeighborTempM1= 95 -- used by MultiTry's "Neighbor Generation" subroutines
QS_TopLevelTempMT= 96 -- used by MultiTry top-level manager
QS_TopLevelTempM2= 97 -- used by MultiTry top-level manager
QS_NeighborTempSA= 98 -- used by SA's "Neighbor Generation" subroutines
QS_TopLevelTempSA= 99 -- used by SA top-level manager

--------------- THINGS THE PROGRAM LEARNS AND REMEMBERS --------------
PuzzleAllowsMutate = false
PuzzleHasLockedSegs = false
PuzzleHasContactMap = false
PuzzleUsesFilters = false
PuzzleIsDesigner = false
ChangedSecondaryStructure = false
ChangedFrozenness = false
EndCalled = false
StartTime = 0
InitialScore = 0.0
ContactMapScore = 0.0
InitialBandCount = band.GetCount( )
CysteineBandPairCount = 0
InitialClashImportance = behavior.GetClashImportance()
HasContactMapSegList = {}
InitialFrozenTable = {}
IsSlowFilterDisabled = behavior.GetSlowFiltersDisabled( )

-- Usage of following tables: bounds = allHelices[i]  use bounds.sIdx and bounds.eIdx
allHelices = { {sIdx = 0, eIdx = 0} }
allSheets = { {sIdx = 0, eIdx = 0} }
allLoops = { {sIdx = 0, eIdx = 0} }
possibleHelices = { {sIdx = 0, eIdx = 0} }
cysteineSegs = {}
cysteinePairs = {}


NonLockedSegList = {}
SecStructList = {}
segCt = structure.GetCount( )
MinUnlockedSeg = 1
MaxUnlockedSeg = segCt
NoMutateList = {}
ShakeAfterMutateList = {}   -- generated from NoMutateList, depending on options wanted by user


---------------- STATISTICS-TRACKING ---------------------
algName = "none"   -- for tracking statistics about how successful an operation is
algNameList = {}
algStats = {
    name = "none",
    posScoreDelta = 0.0, 
    scoreDelta = 0.0, 
    runTime = 0, 
    callCount = 0,
    failCount = 0,
}
TotalShakeTime = 0



----------------------------------------------------------------
-- SIMPLE HELPER FUNCTIONS
----------------------------------------------------------------
function DebugPrint( str )
    if DEBUGRUN then print( str ) end
end

function BoolStr( bval )
    if bval then return "true" end
    return "false"
end

function PadString( strIn, lenWanted )
    strOut = strIn
    while string.len( strOut ) < lenWanted do
        strOut = strOut .. " "
    end
    return strOut
end

function PadNumString( numIn, lenWanted )
    strOut = ""..numIn
    while string.len( strOut ) < lenWanted do
        strOut = strOut .. " "
    end
    return strOut
end

-- a function that parses a string of aa indices and returns them as a sorted array.
function ParseAaString( str, AaArray )
    -- fancy parsing here: Input format is: idx1,idx2-idx3,idx4,idx5,idx6-idx7,...
    for w in string.gmatch( str, "[^,]+" ) do
        if string.find( w, "-" ) ~= nil then
            local s, f = string.match( w, "(%d+)%-(%d+)" )
            start = tonumber( s )
            fin = tonumber( f )
            for i=start,fin do 
                local okToAdd = true
                for j=1,#AaArray do
                    if AaArray[ j ] == i then 
                        okToAdd = false
                        break
                    end
                end
                if okToAdd then 
                    AaArray[ #AaArray+1] = i 
                end
            end
        else
            local seg = tonumber( w )
            local okToAdd = true
            for j=1,#AaArray do
                if AaArray[ j ] == seg then 
                    okToAdd = false
                    break
                end
            end
            if okToAdd then 
                AaArray[ #AaArray+1] = seg 
            end
        end
    end
    if #AaArray > 0 then table.sort(AaArray) end
end

function TrimNum( val )
    return val - val % 0.001 
end

----------------------------------------------------------------
-- ACTION STATISTICS AND PREFERENCES
----------------------------------------------------------------
function WhereIsStringInArray( array, str )
    for i=1, #array do if str == array[ i ] then return i end end
    return 0
end

function AddActionToStats( scoreDelta, timeDelta )
    if WhereIsStringInArray( algNameList, algName ) == 0 then
        algNameList[ #algNameList + 1] = algName
    end

    if algStats[ algName ] == nil then 
        algStats[ algName ] = { 
                  posScoreDelta = 0.0, 
                  scoreDelta = 0.0, 
                  runTime = 0, 
                  callCount = 0,
                  gainCount = 0,
                  failCount = 0,
        }
    end

    if scoreDelta > 0.0 then
        algStats[ algName ].posScoreDelta = algStats[ algName ].posScoreDelta + scoreDelta
        algStats[ algName ].gainCount = algStats[ algName ].gainCount + 1
    end
    algStats[ algName ].scoreDelta = algStats[ algName ].scoreDelta + scoreDelta
    algStats[ algName ].runTime = algStats[ algName ].runTime + timeDelta
    algStats[ algName ].callCount = algStats[ algName ].callCount + 1
end

function AddFailToStats( )
    if WhereIsStringInArray( algNameList, algName ) == 0 then
        algNameList[ #algNameList + 1] = algName
    end

    if algStats[ algName ] == nil then 
        algStats[ algName ] = { 
                  posScoreDelta = 0.0, 
                  scoreDelta = 0.0, 
                  runTime = 0, 
                  callCount = 0,
                  gainCount = 0,
                  failCount = 1,
        }
        return
    end

    algStats[ algName ].failCount = algStats[ algName ].failCount + 1
end

function ReportStats( )
    -- scoreDelta turns out to be not very interesting to look at
    print("Action stats:  " .. #algNameList )
    for i = 1, #algNameList do
        local t = algNameList[ i ]
        if algStats[ t ] ~= nil and algStats[ t ].callCount > 0 then
            print( "  " .. PadString(t, 20) .. ":" .. 
                   " posScore=" .. TrimNum( algStats[ t ].posScoreDelta ) ..
                   " runTime=" .. algStats[ t ].runTime ..
                   " runs=" .. algStats[ t ].callCount ..
                   " gains=" .. algStats[ t ].gainCount ..
                   " fails=" .. algStats[ t ].failCount )
        end
    end
end

function ResetStats( )
    for i = 1, #algNameList do
        local t = algNameList[ i ]
        if algStats[ t ] ~= nil then
            algStats[ t ].scoreDelta = 0
            algStats[ t ].posScoreDelta = 0
            algStats[ t ].runTime = 0
            algStats[ t ].callCount = 0
            algStats[ t ].gainCount = 0
            algStats[ t ].failCount = 0
        end
    end
end

function ChooseAlgorithm( algList )
    -- we dither runs and gains by 1 partly to avoid divide-by-zero, and 
    -- partly to encourage bringing in not-yet-used algorithms.
    local scaleValue = 0.0
    local funcValues = {}

    for i = 1, #algList do
        local runs = 0
        local gains = 0
        if algStats[ algList[ i ] ] == nil then 
            runs = 1
            gains = 1
        else
            -- failCount/2 is a bit artificial, but failed actions should count against an algorithm too.
            runs = 1 + algStats[ algList[ i ] ].callCount + math.floor( algStats[ algList[ i ] ].failCount / 2 )
            gains = 1 + algStats[ algList[ i ] ].gainCount
        end
        funcValues[ i ] = gains / runs    -- finding value of this item
        scaleValue = scaleValue + funcValues[ i ]   -- finding total size of bucket
    end

    local which = random( )
    local sum = 0.0
    for i = 1, #funcValues do
        sum = sum + ( funcValues[ i ] / scaleValue )
        if which < sum then 
            return algList[ i ] 
        end
    end
    return algList[ #algList ]
end

----------------------------------------------------------------
-- BASIC SCORE FUNCTIONS
----------------------------------------------------------------
-- Not for "external" use - call getScore. This could change if customers want 
-- something besides current or recentbest.
function internalGetScore( wantRB, useFilter )
    if wantRB == nil then wantRB = false end
    if useFilter and KILLING_FILTERS then TurnOnSlowFilters( ) end
    local s=0.0
    if not USENORMALSCORE then
        if wantRB then s = recentbest.GetEnergyScore( )
        else s=current.GetEnergyScore( )
        end
    else
        if wantRB then s = recentbest.GetScore( )
        else s=current.GetScore( )
        end
    end
    if useFilter and KILLING_FILTERS then TurnOffSlowFilters( ) end
    return s
end
function getScore( )
    return internalGetScore( false, true )
end
function getQuickScore( )
    return internalGetScore( false, false )
end
function getRBScore( )
    return internalGetScore( true, true )
end

function getSlotScore( qs )
    save.Quicksave( QS_Swapslot )
    save.Quickload( qs )
    sc = getScore( )
    save.Quickload( QS_Swapslot )   
    return sc
end

function SaveBest( forceBest, qs_best )  -- most callers will leave both as nil.
    if qs_best == nil then qs_best = QS_Best end
    if forceBest == nil then forceBest = false end
    
    local score = getScore( )
    local CurrentBestScore = getSlotScore( qs_best )
    
    -- use "score >=" because same-with-unfrozen/unbanded is better
    if forceBest or score >= CurrentBestScore then
        save.Quicksave( qs_best )   
        CurrentBestScore = score
    end

    if PuzzleHasContactMap then
        CurrentBestCmScore = getSlotScore( QS_CmBest )
        cmScore = GetContactScore( )
        if forceBest or (CONTACTMAP_TAKEDROPS and cmScore > CurrentBestCmScore)
        then
            save.Quicksave( QS_CmBest )
            CurrentBestCmScore = score
        end
        if CONTACTMAP_TAKEDROPS and CurrentBestScore ~= CurrentBestCmScore then
            DebugPrint(" In Best("..qs_best.."): score="..CurrentBestScore )
            DebugPrint(" In CmBest("..QS_CmBest.. "): score="..CurrentBestCmScore )
        end
    end
end

function LoadBest( qs_best )
    if qs_best == nil then qs_best = QS_Best end
    local CurrentBestScore = getSlotScore( qs_best )
    
    if PuzzleHasContactMap and CONTACTMAP_TAKEDROPS then
        local CurrentBestCmScore = getSlotScore( QS_CmBest )
        if CurrentBestCmScore >= CurrentBestScore - CONTACTMAP_MAX_DROP
        then
            save.Quickload( QS_CmBest )
        else
            save.Quickload( qs_best )
        end
    else
        save.Quickload( qs_best )
    end
end

function GetBestScore( )
    local CurrentBestScore = getSlotScore( qs_best )
    if PuzzleHasContactMap and CONTACTMAP_ALLOWDROP then
        CurrentBestCmScore = getSlotScore( QS_CmBest )
        if CurrentBestCmScore >= CurrentBestScore - CONTACTMAP_MAX_DROP 
        then
            return CurrentBestCmScore
        else
            return CurrentBestScore
        end
    else
        return CurrentBestScore
    end
end

----------------------------------------------------------------
-- RANDOM FUNCTIONS
----------------------------------------------------------------
function seedRandom( )
    seed=os.time( )/math.abs( getQuickScore( ) )
    seed=seed%0.001
    seed=1/seed
    while seed<10000000 do seed=seed*1000 end
    seed=seed-seed%1
    DebugPrint( "Seed is: "..seed )
    math.randomseed( seed )

    -- throw away a couple of randoms
    math.random( )
    math.random( )
end

function random( n1,n2, forceFloat ) --random function returns int or float depends on input vars
    if forceFloat == nil then forceFloat = false end
    if n1==nil then 
        return math.random( )
    else
        if n2==nil then
            if n1 == 0 then return 0 end -- a random number between 0 and 0 is 0
            if n1%1==0 then -- can't test for "forceFloat", so caller must beware
                return math.random( n1) --integer
            else
                return math.random( ) * n1 --float
            end
        else
            if n1%1==0 and n2%1==0 and not forceFloat then
                return math.random( n1, n2 ) --integer between 
            else
                return math.random( ) * (n2 - n1) + n1 --float between
            end
        end
    end
end

function randomBool( probTrue )
    if probTrue == nil or probTrue == 0.5 then
        local r = random(0, 1, false )
        return r == 0 
    else
        local r = random()
        if r < probTrue then return true 
        else return false 
        end
    end
end

function randomDice( ctDice, minValPerDie, maxValPerDie )  -- distro that prefers "middle" values
    total = 0
    for i=1, ctDice do
        total = total + random( minValPerDie, maxValPerDie )
    end
    return total
end

function randomThetaPhi( )
    return   math.acos( random( -1.0, 1.0, true ) ),  random( 2 * math.pi )
end

function randomizeIndexList( idxList )
    for i=1, #idxList do
        j = random( #idxList )
        if j ~= i then
            idxList[i], idxList[j] = idxList[j], idxList[i]
        end
    end
end

function randomSeg( )
    return random( segCt )
end

function randomUnlockedSeg( )
    if not PuzzleHasLockedSegs then 
        return randomSeg( ) 
    end
    return NonLockedSegList[ random( #NonLockedSegList ) ]
end

function randomAASeg( whichAA )
    local segs = {}
    for i=1, #NonLockedSegList do
        if structure.GetAminoAcid( i ) == whichAA then
            segs[ #segs + 1 ] = i
        end
    end
    if  #segs == 0 then return 0 end
    return segs[ random( #segs ) ]
end

function randomMovableSeg( )
    for i=1, 20 do
        local seg = randomUnlockedSeg( )
        if seg == nil then break end
        local bb, sc = freeze.IsFrozen( seg )
        if not bb then return seg end
    end

    -- if we are here, then we had a hard time finding a nonfrozen seg. 
    -- make an array and search the hard way.
    local movables = {}
    for i=1, #NonLockedSegList do
        local bb, sc = freeze.IsFrozen( seg )
        if not bb then movables[ #movables + 1 ] = i end
    end
    if #movables == 0 then return 1 end -- nothing we can do
    return movables[ random( #movables ) ]
end


function randomMovableFarSeg( segOrigin, minSep, minDist, maxDist )
    local segList = {}
    for e=1, segCt do
        if IsMovableSeg( e ) and (e <= segOrigin - minSep or e >= segOrigin + minSep ) then
            dist = structure.GetDistance( segOrigin, e )
            if dist >= minDist and dist <= maxDist then segList[ #segList + 1 ] = e end
        end
    end
    if #segList >= 1 then return segList[ random( #segList ) ]
    else return 0
    end
end

function randomLowScoringSeg( ctIndicesToChooseAmong, subscore )
    local idxList = {}
    getWorstScoringAAs( idxList, ctIndicesToChooseAmong, subscore )
    return idxList[ random(  #idxList ) ]
end

function randomContactableSeg( )
    return HasContactMapSegList[ random( #HasContactMapSegList ) ]
end

function randomContactSeg( segIn, heatThreshold )
    local segOpts = {}
    for i = 1, segCt do
        if i < segIn - 1 or i > segIn + 1 then
            if contactmap.GetHeat( segIn, i ) >= heatThreshold then
                segOpts[ #segOpts + 1] = i
            end
        end
    end
    if #segOpts == 0 then return 0 end
    return segOpts[ random( #segOpts ) ]
end

function RandomPrimeNotMoreThan( max )
    local primes = {2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97,101}
    lastIdx = #primes
    for i=1, #primes do
        if primes[ i ] > max then break end
        lastIdx = i
    end
    return primes[ math.random( lastIdx ) ]
end

function Coprime( n )
   -- find the highest prime < 70% of "n" that is coprime with "n"
    local primes = {2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97,101}
    for i = #primes, 1, -1 do
        if primes[ i ] < 0.70*n and n % primes[ i ] ~= 0 then
            return primes[ i ]
        end
    end
    return 1
end

----------------------------------------------------------------
-- SCORE-USING FUNCTIONS
----------------------------------------------------------------
function sortAllByScore(segs, scores)
    table.sort(segs, function(n1,n2) return scores[n1] < scores[n2] end)
end

function ConstructScoreMatrix( scmat )
    -- TODO: Get the list of "active" subscores from puzzle info; use that info here.
    for i=1, segCt do
        scmat[i] = {}
        scmat[i][1]   =  current.GetSegmentEnergySubscore( i, "backbone" )
        scmat[i][2]   =  current.GetSegmentEnergySubscore( i, "sidechain" )
        scmat[i][3]   =  current.GetSegmentEnergySubscore( i, "clashing" )
        scmat[i][4]   =  current.GetSegmentEnergySubscore( i, "packing" ) 
        scmat[i][5]   =  current.GetSegmentEnergySubscore( i, "density" )
        scmat[i][6]   =  current.GetSegmentEnergySubscore( i, "hiding" )
        scmat[i][7]   =  current.GetSegmentEnergySubscore( i, "bonding" )
        scmat[i][8]   =  current.GetSegmentEnergySubscore( i, "disulfides" )
        scmat[i][9]   =  current.GetSegmentEnergySubscore( i, "ideality" )
        scmat[i][10]  =  current.GetSegmentEnergySubscore( i, "pairwise" )
        scmat[i][11]  =  current.GetSegmentEnergySubscore( i, "other" )
    end
end

function GetEuclideanDistanceVector( scmat1, scmat2, vecOut )
    score = 0.0
    for j=1, #scmat1[1] do
        vecOut[j] = 0.0
    end
    for i=1, segCt do
        for j=1, #scmat1[i] do
            vecOut[j] = vecOut[j] + (scmat2[i][j] - scmat1[i][j]) * (scmat2[i][j] - scmat1[i][j])
        end
    end
    for j=1, #vecOut do
        vecOut[j] = math.sqrt( vecOut[j] )
    end
end
function IsNeighborDifferentFromOld( oldScoreMat, threshold )
    local newScoreMat = {}
    local scoreVector = {}
    ConstructScoreMatrix( newScoreMat )
    GetEuclideanDistanceVector( oldScoreMat, newScoreMat, scoreVector )
    for i=1, #scoreVector do
        if scoreVector[i] > threshold then return true end
    end
    return false
end

function IsProteinChanged( oldQuickScore, oldScoreMat, threshold )
    if (not KILLING_FILTERS) and (not current.AreConditionsMet( )) then return false end
    local currScore = getQuickScore( )
    if currScore > oldQuickScore then return true end
    return IsNeighborDifferentFromOld( oldScoreMat, threshold )
end

-- returns all aas in order of score (if subscore is nil then does GetSegmentEnergyScore)
function getWorstScoringAAs( idxList, maxWanted, subscore )
    local scoreList = {}
    local idxes = {}
    for i=1, segCt do
       idxes[i] = i
    end
    for i=1, segCt do
        if subscore == nil then
            scoreList[i] = current.GetSegmentEnergyScore(i)
        else
            scoreList[i] = current.GetSegmentEnergySubscore( i, subscore )
        end
    end
    sortAllByScore( idxes, scoreList )
    for i = 1, math.min( maxWanted, segCt ) do
        idxList[i] = idxes[i]
    end
end
function getBestScoringAAs( idxList, maxWanted, subscore )
    local scoreList = {}
    local idxes = {}
    for i=1, segCt do
       idxes[i] = i
    end
    for i=1, segCt do
        if subscore == nil then
            scoreList[i] = current.GetSegmentEnergyScore(i)
        else
            scoreList[i] = current.GetSegmentEnergySubscore( i, subscore )
        end
    end
    sortAllByScore( idxes, scoreList )
    for i = 1, math.min( maxWanted, segCt ) do
        idxList[i] = idxes[segCt - i + 1]
    end
end
function isBadScorer(idx, mustBeAtLeastThisBad)
    local idxList = {}
    getWorstScoringAAs( idxList, mustBeAtLeastThisBad)
    for i=1, #idxList do
        if idxList[i] == idx then return true end
    end
    return false
end

function SortIdxListBySubscore( idxList, subscoreName )
    local subscoreList = {}
    for i=1, #idxList do
        subscoreList[ i ] =  current.GetSegmentEnergySubscore( idxList[ i ], subscoreName )
    end
    sortAllByScore( idxList, subscoreList )
end

function GetWorstScoringRegions( idxList, len )
    scoreList = { }
    for i=1, segCt - len + 1 do
        scoreList[ i ] = 0.0
        idxList[ i ] = i
        for j = 0, len - 1 do
            scoreList[ i ] = scoreList[ i ] + current.GetSegmentEnergyScore( i + j )
        end
    end
    sortAllByScore( idxList, scoreList )
end

----------------------------------------------------------------
-- DISTANCE FUNCTIONS
----------------------------------------------------------------
function computeDistanceSums(distList)
    for i=1, segCt do
        distList[i] = 0.0 
    end
    for i=1, segCt-1 do
        for j=i+1, segCt do
            distList[i] = distList[i] + structure.GetDistance( i, j ) / segCt
            distList[j] = distList[j] + structure.GetDistance( i, j ) / segCt
        end
    end
end
function getCentralAAs( idxList )
    local distList = {}
    for i=1, segCt do
       idxList[i] = i
    end
    computeDistanceSums( distList )
    sortAllByScore( idxList, distList )
end

-- returns all hydrophobics in order of "centralness"
function getCentralHydrophobics(idxList, ctWanted)
    local idxListInternal = {}
    getCentralAAs(idxListInternal)
    for i=1, #idxListInternal do
        if structure.IsHydrophobic( idxListInternal[ i ] ) then
            idxList[ #idxList + 1 ] = idxListInternal[ i ]
            if #idxList == ctWanted then return end
        end
    end
end
function isCentralHydrophobic( idx, mustBeAtLeastThisGood )
    local idxList = {}
    getCentralHydrophobics( idxList, mustBeAtLeastThisGood )
    for i=1, #idxList do
        if idxList[i] == idx then return true end
    end
    return false
end

function getSegsWithinSphere( idxList, centerIdx, radius, includeCenter )
    if includeCenter then
        idxList[ #idxList + 1 ] = centerIdx
    end
    for i=1, segCt do
        if i ~= centerIdx and structure.GetDistance( centerIdx, i ) <= radius then
            idxList[ #idxList + 1 ] = i
        end
    end
end

function getSegsOutsideSphere( idxList, centerIdx, radius )
    for i=1, segCt do
        if i ~= centerIdx and structure.GetDistance( centerIdx, i ) > radius then
            idxList[ #idxList + 1 ] = i
        end
    end
end

function getCMSegsWithinSphere( idxList, centerIdx, radius, includeCenter, heatThreshold )
    if includeCenter then
        idxList[ #idxList + 1 ] = centerIdx
    end
    for i=1, segCt do
        if i ~= centerIdx and structure.GetDistance( centerIdx, i ) <= radius and
            contactmap.GetHeat( centerIdx, i ) >= heatThreshold
        then
            idxList[ #idxList + 1 ] = i
        end
    end
end

function getCMSegsOutsideSphere( idxList, centerIdx, radius, heatThreshold )
    for i=1, segCt do
        if i ~= centerIdx and structure.GetDistance( centerIdx, i ) > radius and
            contactmap.GetHeat( centerIdx, i ) >= heatThreshold
        then
            idxList[ #idxList + 1 ] = i
        end
    end
end

function getSegsWithinRangeOfTube( idxList, startSeg, endSeg, radius )
    -- if a seg is "near" to one of the ones in our [startSeg, endSeg] range, then we want it
    for i=1, segCt do
        for j=startSeg, endSeg do
            if structure.GetDistance( i, j ) < radius then
                idxList[ #idxList + 1] = i
                break
            end
        end  -- for j
    end      -- for i
end

-- produces a list of all AAs "near" the ones in a given inputList
function getSegsNearListed( inputList, radius)
    local retval = {}
    for i=1, segCt do
        for j=1, #inputList do
            if structure.GetDistance(i,j)<radius then retval[ #retval+1 ] = j break end
        end
    end
    return retval
end

function IsLocalMinimum( segCenter, segCheck)
    if segCheck == segCenter then return false end  -- exclude selfies
    ckDist = structure.GetDistance(segCenter, segCheck)
    if segCheck > 1 then 
        ckDist1 = structure.GetDistance( segCenter, segCheck - 1 )
        if ckDist1 < ckDist then return false end
        if segCheck > 2 then 
            ckDist2 = structure.GetDistance( segCenter, segCheck - 2 )
            if ckDist2 < ckDist then return false end
        end
    end
    if segCheck < segCt then
        ckDist1 = structure.GetDistance( segCenter, segCheck + 1 )
        if ckDist1 < ckDist then return false end
        if segCheck < segCt - 1 then 
            ckDist2 = structure.GetDistance( segCenter, segCheck + 2 )
            if ckDist2 < ckDist then return false end
        end
    end
    return true
end

function IsLocalMaximum( segCenter,  segCheck)
    if segCheck == segCenter then return false end  -- exclude selfies
    ckDist = structure.GetDistance(segCenter, segCheck)
    if segCheck > 1 then 
        ckDist1 = structure.GetDistance( segCenter, segCheck - 1 )
        if ckDist1 > ckDist then return false end
        if segCheck > 2 then 
            ckDist2 = structure.GetDistance( segCenter, segCheck - 2 )
            if ckDist2 > ckDist then return false end
        end
    end
    if segCheck < segCt then
        ckDist1 = structure.GetDistance( segCenter, segCheck + 1 )
        if ckDist1 > ckDist then return false end
        if segCheck < segCt - 1 then 
            ckDist2 = structure.GetDistance( segCenter, segCheck + 2 )
            if ckDist2 > ckDist then return false end
        end
    end
    return true
end

function GetNearestPhobicsOutsideRange( seg, segExcludeStart, segExcludeEnd, nearPhobics, ctWanted )
    local idxListInternal = {}
    local distList = {}
    for i=1, segCt do idxListInternal[ #idxListInternal + 1 ] = i end

    for i=1, segExcludeStart - 1 do
         if structure.IsHydrophobic( i ) then 
             distList[ #distList + 1 ] = structure.GetDistance( seg, i )
         else
             distList[ #distList + 1 ] = 1000
         end
    end
    for i=segExcludeStart, segExcludeEnd do
        distList[ #distList + 1 ] = 10000
    end
    for i=segExcludeEnd + 1, segCt do
         if structure.IsHydrophobic( i ) then 
             distList[ #distList + 1 ] = structure.GetDistance( seg, i )
         else
             distList[ #distList + 1 ] = 1000
         end
    end

    sortAllByScore( idxListInternal, distList )

    for i=1, ctWanted do
        if i > #idxListInternal then break end
        nearPhobics[ #nearPhobics + 1 ] = idxListInternal[ i ]
    end
end

function getSegsNotInList( list )
    table.sort( list )
    local result={ }
    if #list > 0 then for j = 1, list[1]-1 do result[ #result+1 ] = j end end
    for i = 1, #list-1 do
       for j = list[ i ]+1, list[ i+1 ]-1 do result[ #result+1 ]=j end
    end
    for j = list[ #list ]+1, segCt do result[ #result+1 ]=j end
    return result
end

----------------------------------------------------------------
-- SIDECHAIN AND ATOM FUNCTIONS
----------------------------------------------------------------
aaCount = 20
aaPhilCount = 9
aaPhobCount = 11
function AA(idx)
    -- hydrophilics
    if idx == 1 then return "r"
    elseif idx == 2 then return "s"
    elseif idx == 3 then return "t"
    elseif idx == 4 then return "n"
    elseif idx == 5 then return "d"
    elseif idx == 6 then return "q"
    elseif idx == 7 then return "e"
    elseif idx == 8 then return "h"
    elseif idx == 9 then return "k"
    -- hydrophobics
    elseif idx == 10 then return "g"
    elseif idx == 11 then return "a"
    elseif idx == 12 then return "c"
    elseif idx == 13 then return "v"
    elseif idx == 14 then return "l"
    elseif idx == 15 then return "i"
    elseif idx == 16 then return "m"
    elseif idx == 17 then return "p"
    elseif idx == 18 then return "f"
    elseif idx == 19 then return "y"
    else return "w" end      -- idx == 20
end

function AAphilic(idx)
    -- hydrophilics
    if idx == 1 then return "r"
    elseif idx == 2 then return "s"
    elseif idx == 3 then return "t"
    elseif idx == 4 then return "n"
    elseif idx == 5 then return "d"
    elseif idx == 6 then return "q"
    elseif idx == 7 then return "e"
    elseif idx == 8 then return "h"
    else return "k"      -- idx == 9
    end
end

function AAphobic(idx)
    -- hydrophobics
    if idx == 1 then return "g"
    elseif idx == 2 then return "a"
    elseif idx == 3 then return "c"
    elseif idx == 4 then return "v"
    elseif idx == 5 then return "l"
    elseif idx == 6 then return "i"
    elseif idx == 7 then return "m"
    elseif idx == 8 then return "p"
    elseif idx == 9 then return "f"
    elseif idx == 10 then return "y"
    else return "w"      -- idx == 11
    end 
end

-- for branched aas, simply picks one (longer, one with donor/acceptor tip, or if no difference then either)
function GetAtomOfTip( aa )
    if aa == "a" then return 10
    elseif aa == "c" then return 10 -- real tip - 1 because maybe cysteine bridge
    elseif aa == "d" then return 8
    elseif aa == "e" then return 9
    elseif aa == "f" then return 20
    elseif aa == "g" then return 0 -- glycine has no tip. just use a backbone atom
    elseif aa == "h" then return 17
    elseif aa == "i" then return 18
    elseif aa == "k" then return 9
    elseif aa == "l" then return 16
    elseif aa == "m" then return 17
    elseif aa == "n" then return 14
    elseif aa == "p" then return 13
    elseif aa == "q" then return 9
    elseif aa == "r" then return 22
    elseif aa == "s" then return 11
    elseif aa == "t" then return 6
    elseif aa == "v" then return 13
    elseif aa == "w" then return 24
    elseif aa == "y" then return 12
    else return 0
    end
end

function GetTipAtomOfSeg( idx )
	-- cysteine loses last atom if disulfide bridge
	if structure.GetAminoAcid( idx ) == "c" then
		return structure.GetAtomCount( idx )
	else
		return GetAtomOfTip( structure.GetAminoAcid( idx ))
	end
end

function GetCysteineSulfurAtom( idx )
	if structure.GetAminoAcid(idx) ~= "c" then
		return 0 
	end
	totalAtoms = structure.GetAtomCount(idx)
	if totalAtoms == 10 then 
		-- there is a cysteine bond, so things shifted
		return 7
	else
		return 6
	end
end

function PickAtomNumberForBand( idx )
    if not BAND_TO_SIDECHAINS then return CENTER_CARBON end

    local r = random( 1, 7 )   -- consider adjusting probability?
    if r <= 3 then         -- do this for r = 1,2,3
        return BETA_CARBON
    elseif r == 4 then   -- do this rarely, because it's only rarely useful
        return GetTipAtomOfSeg( idx )
    else                      -- do this for r = 5,6,7
        return CENTER_CARBON
    end
end

function SetCysteineSegList( )
    cysteineSegs = {}
    for i=1, segCt do 
        lt = structure.GetAminoAcid( i )
        if lt == "c" or lt == "C" then 
            cysteineSegs[#cysteineSegs + 1] = i 
        end
    end
end

function GuessCysteinePairs( )
    ctPairsWanted = math.floor( #cysteineSegs / 2 )
    cysteinePairs = {}
    localCystList = {}
    for i=1, #cysteineSegs do
        localCystList[i] = cysteineSegs[i]
    end
    for ii = 1, ctPairsWanted do
        minDist = 1000000.0
        bestI = 0
        bestJ = 0
        for i = 1, #localCystList - 1 do
            cyst1 = localCystList[i]
            for j=i+1, #localCystList do
                cyst2 = localCystList[j]
                -- don't try to match up "adjacent" cysteines.
                if cyst2 - cyst1 > 3 then 
                    distBetween = structure.GetDistance( cyst1, cyst2 )
                    if distBetween < minDist then
                        bestI = i
                        bestJ = j
                        minDist = distBetween
                    end
                end
            end
        end
        cysteinePairs[#cysteinePairs+1] = {localCystList[bestI], localCystList[bestJ]}
        table.remove( localCystList, bestJ )
        table.remove( localCystList, bestI )
    end
    return true
end

function DoesCysteineBondExist( idx1, idx2)		
	-- NOTE: I have problems if two cysteine bonds exist involving the
	--       two given segments but the bonds are not actually between 
	--       my two segments (eg a<->c and b<->d exist, I want a<->b).
	--
	--       My attempts to use distance are kludgy but I don't have
	--       better (can't ask for distance between atoms), so I will 
	--       leave the problem in situ.
	--       structure.GetAtomCount could also be checked (10 vs 11)
    --       but that doesn't help the core problem here.
	
	-- dist <= 6.996 seems to correspond to success.
	dist = structure.GetDistance(idx1, idx2 )
	dScore1 = current.GetSegmentEnergySubscore( idx1, "disulfides" )
	dScore2 = current.GetSegmentEnergySubscore( idx2, "disulfides" )
	bondExists = (dist < 7.5 and dScore1 ~= 0.0 and dScore2 ~= 0.0)
	
    --if bondExists then
    --    print( "CBond "..idx1.."-"..idx2..
	--		   " dist="..TrimNum(dist)..
	--		   "   sc1="..TrimNum(dScore1)..
	--		   "   sc2="..TrimNum(dScore2))
	--end
	
	return bondExists
end

function CheckCysteineSanity()
	for i = 1, #cysteinePairs do
		pair = cysteinePairs[i]
		dist = structure.GetDistance( pair[1], pair[2] )
		if dist > MAX_RI_CYSTEINE_SEPARATION then
			print( "Rejected action - cysteine "..pair[1].." far from "..pair[2])
			return false
		end
		--print("distance "..pair[1].." to "..pair[2].." is "..dist)
	end
	return true
end


----------------------------------------------------------------
-- MUTATION-RELATED FUNCTIONS
----------------------------------------------------------------
function checkAla( )
    local ala=0
    for n = 1, segCt do
        segType = structure.GetAminoAcid(n)
        if segType == 'a' or segType == 'g' then
            ala=ala+1
        end
    end
    return ala == segCt
end

function IsPuzzleMutable( )
   for i=1, segCt do
      if structure.IsMutable(i) then return true end
   end
   return false
end

function IsNoMutate( idx )
   if not ALLOW_SPOT_MUTATION then return true end
   if not PuzzleAllowsMutate then return true end
   if not structure.IsMutable(idx) then return true end

   for i=1,#NoMutateList do
       if idx == NoMutateList[i] then return true end
   end
   return false
end

----------------------------------------------------------------
-- SLOW FILTER FUNCTIONS
----------------------------------------------------------------
function TurnOnSlowFilters( )
    if behavior.GetSlowFiltersDisabled( ) then
        --print("slow filters are now ON")
        behavior.SetSlowFiltersDisabled( false )
    end
end

function TurnOffSlowFilters( )
    if not behavior.GetSlowFiltersDisabled( ) then
        --print("slow filters are now OFF")
        behavior.SetSlowFiltersDisabled( true )
    end
end

function ResetSlowFilterState( )
    --if IsSlowFilterDisabled then 
    --    print("slow filters are now ON")
    --end
    behavior.SetSlowFiltersDisabled( IsSlowFilterDisabled )
end

function CheckForSlowFilters( )
    TurnOffSlowFilters( )
    local retval = not current.AreConditionsMet()
    ResetSlowFilterState( )
    return retval
end

----------------------------------------------------------------
-- LOCKED/FROZEN FUNCTIONS
----------------------------------------------------------------
function GenerateNonLockedSegList( )
    for i=1, segCt do
        if structure.IsLocked( i ) then 
            PuzzleHasLockedSegs = true
        else
            NonLockedSegList[ #NonLockedSegList + 1] = i
        end
    end
end

function IsUnlockedSeg( seg )
   return not structure.IsLocked( seg )
end

function FindAllMovableSegs( )
    if PuzzleHasLockedSegs then
        for i=1, segCt do
            if IsUnlockedSeg( i ) then 
                MinUnlockedSeg = i
                break
            end
        end
        for i=segCt, 1, -1 do
            if IsUnlockedSeg( i ) then 
                MaxUnlockedSeg = i
                break
            end
        end
    end

    return false
end

function IsMovableSeg( seg1 )
   if seg1 == nil or seg1 < 1 or seg1 > segCt then return false end

   local BB, SC = freeze.IsFrozen( seg1 )
   local sl = structure.IsLocked( seg1 )
   return not (BB or SC or sl)
end

function IsAllMovableRange(seg1, seg2)
    if seg1 == nil or seg2 == nil then return false end
    if seg1 > seg2 then seg1, seg2 = seg2, seg1 end
    if seg1 < 1 or seg2 > segCt then return false end
    for i=seg1, seg2 do
        if not IsMovableSeg( i ) then return false end
    end
    return true
end

function SelectAllMovable( )
    for i=1, segCt do
        if IsMovableSeg( i ) then selection.Select( i ) end
    end
end

function SelectNonFrozens( idxList )
    local segs = { }
    local segsToCheck = idxList

    if segsToCheck == nil then
        segsToCheck = { }
        for i = 1, segCt do
            segsToCheck[ #segsToCheck + 1] = i
        end
    end
    selection.DeselectAll( )
    for i = 1, #segsToCheck do
         b,s = freeze.IsFrozen( segsToCheck[ i ] )
         if not s then 
             segs[ #segs + 1 ] = segsToCheck[ i ]
        end
    end
    for i = 1, #segs do
        selection.Select( segs[ i ] )
    end
end

function ResetFrozenness( )
    freeze.UnfreezeAll( )
    for i=1, segCt do
        bFroz, sFroz = InitialFrozenTable[ i ]
        if bFroz or sFroz then
            freeze.Freeze( i, bFroz, sFroz )
        end
    end
    ChangedFrozenness = false
end

function SaveFrozenness( )
    for i=1, segCt do     -- save "frozenness"
        local bFroz, sFroz = freeze.IsFrozen( i )
        InitialFrozenTable[ i ] = bFroz, sFroz
    end
end

----------------------------------------------------------------
-- CONTACT MAP FUNCTIONS
----------------------------------------------------------------
function CheckForContactMap( )
    PuzzleHasContactMap = false
    local saveval = 0.0
    for i=1, segCt-1 do
        for j=i+1, segCt do
            val = contactmap.GetHeat( i, j )
            if saveval  ~= 0.0 and val ~= saveval then
                PuzzleHasContactMap = true
                ContactMapScore = GetContactScore( )
                return   -- all we wanted to know was whether scores exist
            end
            if saveval == 0.0 then saveval = val end
        end
    end
    return
end

function SegHasContactData( segIn, heatThreshold )
    for i = 1, segCt do
        if i < segIn - 1 or i > segIn + 1 then
            if contactmap.GetHeat( segIn, i ) >= heatThreshold then return true end
        end
    end
    return false
end

function InitializeContactMapSegList( heatThreshold )
    HasContactMapSegList = {}
    for i = 1, segCt do
         if SegHasContactData( i, heatThreshold ) then 
             HasContactMapSegList[ #HasContactMapSegList + 1 ] = i 
        end
    end
end

function GetContactScore( )
    if not PuzzleHasContactMap then return 0 end
    local sum = 0.0
    local segCt = structure.GetCount( )
    for i = 1,segCt-1 do
        for j = i + 1, segCt do
            if contactmap.IsContact( i, j ) then sum = sum + contactmap.GetHeat( i, j ) end
         end
    end
    return sum
end

function ActionDamagedContactMap( )
    if PuzzleHasContactMap then 
        local sc = GetContactScore( )
        if sc < ContactMapScore then return true end
        ContactMapScore = sc   -- now we have a "better" score, so keep that one
    end
    return false
end

----------------------------------------------------------------
-- SECONDARY STRUCTURES
----------------------------------------------------------------

function isExcellentHelixSeg( seg )
    aa = structure.GetAminoAcid( seg )
    return aa == 'a' or aa == 'e' or aa == 'm' or aa == 'k' 
end

function isGoodHelixSeg( seg )
    aa = structure.GetAminoAcid( seg )
    return isExcellentHelixSeg(seg) or 
        aa == 'f' or aa == 'j' or aa == 'l' or aa == 'q' or aa == 'v' or aa == 'w'
end

function isOkHelixSeg( seg )
    aa = structure.GetAminoAcid( seg )
    return isGoodHelixSeg(seg) or 
        aa == 'h' or aa == 'r' or aa == 's' or aa == 't'
end

function isHelixBreaker( seg )
    return structure.GetAminoAcid( seg ) == 'p'
end

function isBadNonBreakerHelixSeg( seg )
    aa = structure.GetAminoAcid( seg )
    return aa == 'g' or aa == 'd' or aa == 'n'
end

function LocateHelices( )
    inside = false
    helixStart = 0
    allHelices = {}
    for i=1, segCt do
        if IsUnlockedSeg( i ) then lastUnlocked = i end
        if ( IsUnlockedSeg( i ) and structure.GetSecondaryStructure(i) == "H" ) then
            if ( inside == false ) then
                inside = true
                helixStart = i
            end
        elseif inside == true then
            inside = false
            allHelices[#allHelices+1] = { sIdx = helixStart, eIdx = i-1 }
        end 
    end 

    -- deal with "last seg is helix"
    if ( inside == true ) then
        allHelices[#allHelices+1] = {sIdx = helixStart, eIdx = lastUnlocked }
    end
end

function LocateSheets( )
    inside = false
    sheetStart = 0
    allSheets = {}
    for i=1, segCt do
        if IsUnlockedSeg( i ) then lastUnlocked = i end
        if ( IsUnlockedSeg( i ) and structure.GetSecondaryStructure(i) == "E" ) then
            if ( inside == false ) then
                inside = true
                sheetStart = i
            end
        elseif inside == true then
            inside = false
            allSheets[#allSheets + 1] = {sIdx = sheetStart, eIdx = i - 1}
        end -- if/else 'E'
    end -- for (segCt)

    -- deal with "last seg is sheet"
    if ( inside == true ) then
        allSheets[#allSheets + 1] = {sIdx = sheetStart, eIdx = lastUnlocked}
    end
end

function LocateLoops( )
    inside = false
    loopStart = 0
    allLoops = {}
    for i=1, segCt do
        if IsUnlockedSeg( i ) then lastUnlocked = i end
        if ( IsUnlockedSeg( i ) and structure.GetSecondaryStructure(i) == "L" ) then
            if ( inside == false ) then
                inside = true
                loopStart = i
            end
        elseif inside == true then
            inside = false
            allLoops[ #allLoops + 1 ] = {sIdx = loopStart, eIdx = i - 1}
        end -- if/else 'L'
    end -- for (segCt)

    -- deal with "last seg is loop"
    if ( inside == true ) then
        allLoops[ #allLoops + 1 ] = {sIdx = loopStart, eIdx = lastUnlocked}
    end
end


function FindFirstProbableHelixInSeglist( sIdx, eIdx )
    startHelix = 0
    endHelix = 0

    if eIdx - sIdx < 6 then  -- not worth looking for one here
        return 0, 0 
    end
    for i = sIdx, eIdx - 5 do
        startHelix = 0
        endHelix = 0
        if isGoodHelixSeg( i ) then
            --print("Possible start at idx="..i.." aa="..structure.GetAminoAcid(i) )
            startHelix = i
            for j=i+1, eIdx do
                endHelix = j
                if isOkHelixSeg( j ) then
                    -- reset mediocre count because we don't mind
                    -- isolated mediocres, just don't want 2 in a row
                    --print("Saw good "..structure.GetAminoAcid(j).." at "..j)                  
                elseif isHelixBreaker( j ) then
                    endHelix = j-1
                    --print("Saw proline at "..j)
                    break
                elseif isBadNonBreakerHelixSeg( j ) then
                    endHelix = j-1
                    --print("Saw bad "..structure.GetAminoAcid(j).." at "..j)
                    break
                else 
                    --print( "Saw mediocre "..structure.GetAminoAcid(j).." at "..j )
                    nextOk = true
                    nextNextOk = true
                    if j < segCt - 1 then
                        nextNextOk = isOkHelixSeg( j + 2 )
                    end
                    if j < segCt then 
                        nextOk = isOkHelixSeg( j + 1 )
                    end
                    if not (nextOk and nextNextOk) then
                        break
                    end
                end
            end
            if endHelix - startHelix > 4 then 
                return startHelix, endHelix 
            end
        end
    end

    return 0, 0
end

function LocateMissingHelices()
    possibleHelices = {}
    for i = 1, #allLoops do
        bounds = allLoops[i]
        --print( "looking at loop: "..bounds.sIdx.."-"..bounds.eIdx )
        if bounds.eIdx - bounds.sIdx >= 7 then
            currStart = bounds.sIdx
            while currStart < bounds.eIdx - 4 do
                startHelix, endHelix = FindFirstProbableHelixInSeglist( currStart, bounds.eIdx )
                if startHelix == 0 then
                    break -- no more helices to find in this stretch of loopiness
                else
                    --print("Found possible helix at s,e="..startHelix..","..endHelix)
                    possibleHelices[#possibleHelices + 1] = {sIdx = startHelix, eIdx = endHelix}
                    currStart = endHelix + 2
                end
            end
        end
    end
end

function isPossibleHelix( idx )
    for i = 1, #possibleHelices do
        bound = possibleHelices[ i ]
        if bound.sIdx <= idx and bound.eIdx >= idx then
            s = bound.sIdx
            e = bound.eIdx
            return true
        end
    end
    return false
end

function GetRegionStartAndEnd( idx ) 
    s = idx
    e = idx
    if not IsUnlockedSeg( idx ) then return idx,idx end

    local typeIdx = structure.GetSecondaryStructure( idx )
    if idx > MinUnlockedSeg then
        s = idx - 1
        while s >= MinUnlockedSeg do
            if structure.GetSecondaryStructure( s ) ~= typeIdx or not IsUnlockedSeg( s ) then break end
            s = s - 1
        end
        s = s + 1
    end
    if idx < MaxUnlockedSeg then
        e = idx + 1
        while e <= MaxUnlockedSeg  do
            if structure.GetSecondaryStructure( e ) ~= typeIdx  or not IsUnlockedSeg( e ) then break end
            e = e + 1
        end
        e = e - 1
    end
    
    if typeIdx == 'l' then
        -- if we think there's a helix hidden here, let's try to offer its boundaries
        for i = 1, #possibleHelices do
            bound = possibleHelices[ i ]
            if bound.sIdx <= idx and bound.eIdx >= idx then
                s = bound.sIdx
                e = bound.eIdx
                return s, e
            end
        end
    end
    return s, e
end

-- uses GetRegionStartAndEnd, but tries to tune "loop" regions a bit
function GetRegionForOperation( idx )
    -- try to find an acceptable region (really bizarre choices of SS can give trouble)
    if not IsUnlockedSeg( idx ) then 
        return idx, idx 
    end
    sIdx, eIdx = GetRegionStartAndEnd( idx )

    -- Check: if our entire region is only 1 or 2 segs long, then see about expanding it
    local regionType = structure.GetSecondaryStructure( idx )
    if eIdx - sIdx <= 1 then 
        -- small range. Try to expand it...
        local sIdx2, eIdx2
        if sIdx > 2 then
            if structure.GetSecondaryStructure(sIdx - 1) == regionType then
                sIdx2, eIdx2 = GetRegionStartAndEnd( sIdx - 1)
                sIdx = sIdx2
            end
         end
        if eIdx < segCt - 1 then
            if structure.GetSecondaryStructure( eIdx + 1 ) == regionType then
                sIdx2, eIdx2 = GetRegionStartAndEnd( eIdx + 1 )
                eIdx = eIdx2
            end
        end 
        
        if eIdx - sIdx <= 1 then 
            -- last try: unless we're at an end, just give ourselves one before and after
            if sIdx > 1 then sIdx = sIdx - 1 end
            if eIdx < segCt then eIdx = eIdx + 1 end
        end
        
        if eIdx - sIdx <= 1 then 
            return idx,idx 
        end  -- failed, something is very wrong
    end

    -- Check: if loop and more than 30-long, then only take part of it
    if regionType == 'L' and eIdx - sIdx > 30 and not isPossibleHelix( sIdx ) then
        -- trim it; it is too long
        if idx - sIdx < 15 then
            eIdx = idx + 15  -- trim the other, longer side
        elseif eIdx - idx < 15 then
            sIdx = idx - 15  -- trim the other, longer side
        else
            -- more than 15 segs slop on both sides!
            -- try to be random and "medium-long": 3d5 on each side should do
            sIdx = idx - randomDice(3, 1, 5)
            eIdx = idx + randomDice(3, 1, 5)
        end
    end

    return sIdx, eIdx
end

function ResetSecondaryStructures( min, max )
    selection.DeselectAll( )
    for i=min, max do 
        if SecStructList[ i ] == "l" then
            selection.Select( i )
        end
    end
    structure.SetSecondaryStructureSelected( "l" )

    selection.DeselectAll( )
    for i=min, max do 
        if SecStructList[ i ] == "e" then
            selection.Select( i )
        end
    end
    structure.SetSecondaryStructureSelected( "e" )

    selection.DeselectAll( )
    for i=min, max do 
        if SecStructList[ i ] == "h" then
            selection.Select( i )
        end
    end
    structure.SetSecondaryStructureSelected( "h" )
    selection.DeselectAll( )
end

function StoreSecondaryStructure( )
    for i=1, segCt do 
        SecStructList[ i ] = structure.GetSecondaryStructure( i ) 
    end
end

function SetSS( min, max, type )
    selection.DeselectAll(  )
    selection.SelectRange( min, max )
    structure.SetSecondaryStructureSelected( type ) 
    selection.DeselectAll(  )
end

function SetAllLoop( min, max )
    SetSS( min, max, "l" )
end

----------------------------------------------------------------
-- CLASH IMPORTANCE
----------------------------------------------------------------
function SetCI( ci )
    behavior.SetClashImportance( SCALE_MAXCI * ci )
end

function ExactSetCI( ci )
    behavior.SetClashImportance( ci )
end

function GetCI( )
    return behavior.GetClashImportance( )
end

----------------------------------------------------------------------
-- CUTTING FUNCTIONS 
----------------------------------------------------------------------
function DeleteCutAtIndex( idx )
    if idx == nil or idx < 1 or idx > segCt then return true end
    local sc = getQuickScore( )
    structure.DeleteCut( idx )
    return getQuickScore( ) ~= sc  -- if no score change at all, then cut deletion didn't happen!
end

function CutEverywhereInRange( min, max )
    if min == nil or min > segCt then min = 1 end
    if max == nil or max < min then max = segCt end
    for i=min, max do
        structure.InsertCut( i )
    end
end

function DeleteAllCuts( min, max )
    if min == nil or min > segCt then min = 1 end
    if max == nil or max < min then max = segCt end
    for i=min, max do
        structure.DeleteCut( i )
    end
end

function AddRandomCuts( cutList, maxCuts )
    for i=1, maxCuts do
        local seg = randomMovableSeg( )
        local gotMatch = false
        for j=1, #cutList do 
            if cutList[ j ] == seg then gotMatch = true end 
        end
        if not gotMatch then
            structure.InsertCut( seg )
            cutList[ #cutList + 1 ] = seg
        end
    end
end

function CleanCutsInList( cutList )
    changeSucceeded = true
    for i=1, #cutList do
        changeSucceeded = DeleteCutAtIndex( cutList[ i ] ) and changeSucceeded
    end
    return changeSucceeded
end

function ResetCuts( )
    if not PRESERVE_USER_CUTS then
        for i=1, segCt do DeleteCutAtIndex( i ) end   -- any cuts that appeared better go now, too
    end
end


------------------------------------------------------------------------------
--
---- WIGGLER FUNCTIONS
--
------------------------------------------------------------------------------

------------------------------------------------------------------------------
-- GENERAL-PURPOSE SHAKERS AND WIGGLERS
------------------------------------------------------------------------------
-- no recentbest here because sometimes we'd call while decreasing score
function ShakeOrMutate( iters, idxList )
    local startTime = os.time( )
    if iters == nil then iters = SHAKESIDECHAIN_ITERS end

    if MUTATE_NOT_SHAKE then
       SelectNonFrozens( idxList )
       if #NoMutateList > 0 then
            -- deselect all no-mutate aas
            for i = 1, #NoMutateList do
                selection.Deselect( NoMutateList[ i ] )
            end
        end

        structure.MutateSidechainsSelected( iters )

        if #ShakeAfterMutateList > 0 then
            selection.DeselectAll( )
            for i = 1, #ShakeAfterMutateList do
                local b,s = freeze.IsFrozen( ShakeAfterMutateList[ i ] )
                if not s then selection.Select( ShakeAfterMutateList[ i ] ) end
            end
            structure.ShakeSidechainsSelected( iters )
        end
    else
        SelectNonFrozens( idxList )
        structure.ShakeSidechainsSelected( iters )
    end

    selection.DeselectAll( )
    TotalShakeTime = TotalShakeTime + ( os.time( ) - startTime )
end

function LongWiggle( )
    if USEITERATIVELONGWIGGLE then
        IterativeWiggleAll( true, true )
    else
        structure.WiggleAll( LONGWIGGLE_ITERS )
    end
end

-- no recentbest here because sometimes we'd call while decreasing score
function WiggleRange( startSeg, endSeg, addFreeze, doGlobal )
    local firstSeg = math.max( MinUnlockedSeg, startSeg )
    local lastSeg = math.min( endSeg, MaxUnlockedSeg )
    if lastSeg < firstSeg then return end -- nothing to do

    -- set up for the local wiggle
    selection.DeselectAll( )
    if addFreeze then
        ChangedFrozenness = true
        if firstSeg > MinUnlockedSeg then freeze.Freeze( firstSeg - 1, true, true ) end
        if lastSeg < MaxUnlockedSeg then freeze.Freeze( lastSeg + 1, true, true ) end
    end

    -- actually perform the local wiggle
    selection.SelectRange( firstSeg, lastSeg )
    if doGlobal then
        structure.WiggleSelected( QUICKWIGGLE_ITERS, true, true )
    else
        structure.LocalWiggleSelected( QUICKWIGGLE_ITERS, true, true )
    end

    -- clean up what we did
    ResetFrozenness( )
    selection.DeselectAll( )
end

-- no recentbest here because sometimes we'd call while decreasing score
function IterativeWiggleAll( doBackbone, doSidechain )
    local lastScore = getQuickScore( )
    local gain = QUICKWIGGLE_GAIN
    recentbest.Save( )
    while gain >= QUICKWIGGLE_GAIN do
        structure.WiggleAll( QUICKWIGGLE_ITERS, doBackbone, doSidechain )  -- do's ok for nil
        recentbest.Restore( )
        gain = getQuickScore( ) - lastScore
        lastScore = getQuickScore( )
    end
end

---------------------------------------------
-- SIMPLE STABILIZER WIGGLES
---------------------------------------------
-- currently unused
function IterativeWiggleList( idxList, minGain, doBackbone, doSidechain )
    print( "IterativeWiggleList" )
    selection.DeselectAll( )
    for i=1, #idxList do selection.Select( i ) end

    local lastScore = getQuickScore( )
    local gain = minGain
    while gain >= minGain do
        structure.WiggleSelected( QUICKWIGGLE_GAIN, doBackbone, doSidechain )  -- nil is allowed for do*
        gain = getQuickScore( ) - lastScore
        lastScore = getQuickScore( )
    end
    selection.DeselectAll( )
end

function WiggleShakeWiggleList( idxList )
    print( "WiggleShakeWiggleList" )
    recentbest.Save( )
    sc = getQuickScore( )

    selection.DeselectAll( )
    for i=1, #idxList do selection.Select( i ) end
    structure.WiggleSelected( LONGWIGGLE_ITERS )
    selection.DeselectAll( )
    recentbest.Restore( )

    ShakeOrMutate( SHAKESIDECHAIN_ITERS, idxList )
    recentbest.Restore( )

    if math.abs( getQuickScore( ) - sc) > 0.1 then -- otherwise, why bother?
        LongWiggle( )  
    end
    recentbest.Restore( )
end

function qStabWiggle( doExtraWork )
    print( "qStabWiggle" )
    if doExtraWork == nil then doExtraWork = false end
    SetCI( 1.0 )
    recentbest.Save( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    recentbest.Restore( )

    if doExtraWork then
        SetCI( 0.4 )
        LongWiggle( )  
        recentbest.Restore( )
        SetCI( 1.0)
        ShakeOrMutate( SHAKESIDECHAIN_ITERS )
        recentbest.Restore( )
    end
    SetCI( 1.0)
    LongWiggle( )  
    recentbest.Restore( )
end

-- the only user for now is the idealizer
-- tends to produce transients
function WiggleZlbNicely( strength, skipLower, skipUpper )
    print( "WiggleZlbNicely" )
    if skipLower == nil then skipLower = MaxUnlockedSeg+1 end
    if skipUpper == nil then skipUpper = MinUnlockedSeg-1 end

    recentbest.Save( )
    ZLB(strength, skipLower, skipUpper)
    LongWiggle( )  
    recentbest.Restore( )
    DelBands( )

    recentbest.Save( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    recentbest.Restore( )
end

function SimpleRegionWiggle( idxLo, idxHi, pickGlobalWiggle, wiggleExpansion )
    if pickGlobalWiggle then
        print( "SimpleRegionWiggle-global")
    else
        print( "SimpleRegionWiggle-local")
    end
    if wiggleExpansion == nil then
        wiggleExpansion = 1 + random( math.floor( segCt / 6 ))
    end
    local startIdx = math.max( MinUnlockedSeg, idxLo - wiggleExpansion )
    local endIdx = math.min( MaxUnlockedSeg, idxHi + wiggleExpansion )
    WiggleRange( startIdx, endIdx, false, pickGlobalWiggle )
end

function VerySimpleShakeWiggle( )
    print( "VerySimpleShakeWiggle" )
    recentbest.Save( )
    ExactSetCI( 0.05 )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    SetCI( 1.0 )
    LongWiggle( )  
    recentbest.Restore( )
end

function SimpleWiggleShake( )
    print( "SimpleWiggleShake" )
    recentbest.Save( )
    local oldRbScore = getRBScore( )
    SetCI( COREPULLING_CI )
    if randomBool( PROB_LOCALWIGGLE ) then
        structure.LocalWiggleAll( QUICKWIGGLE_ITERS )
    else
        structure.WiggleAll( QUICKWIGGLE_ITERS )
    end
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    local newRbScore = getRBScore( )
    if newRbScore > oldRbScore then recentbest.Restore( ) end
end

---------------------------------------------------------
--BAND-CHANGING WRAPUP WIGGLERS
---------------------------------------------------------
function CleanShakeQStab()
    DelBands( )
    ResetFrozenness( )    
    SetCI( 1.0 )
    
    -- shake first to try to make modification acceptable
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    -- wiggle to complete the action
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA ))
end

function ShakeWiggleCleanShakeWiggle( )
    print( "ShakeWiggleCleanShakeWiggle" )
    recentbest.Save( )
    local oldRbScore = getRBScore( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    if randomBool( PROB_LOCALWIGGLE) then
        structure.LocalWiggleAll( 1 )
    else
        structure.WiggleAll( 1 )
    end
    local newRbScore = getRBScore( )
    if newRbScore > oldRbScore then recentbest.Restore( ) end
    DelBands( )
    ResetFrozenness( )

    recentbest.Save( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    LongWiggle( )  
    recentbest.Restore( )
end


function BandedWiggleUntilEnoughChange( scoreThreshold, pullingCI )
    print( "BandedWiggleUntilEnoughChange" )
    if pullingCI == nil then pullingCI = COREPULLING_CI end
    
    SetCI( pullingCI )
    local ss=getQuickScore( )
    for str=0.30, 1.50, 0.11 do   --search enough band strength to move
        SetMyBandStrengths( str )
        if randomBool( PROB_LOCALWIGGLE ) then
            structure.LocalWiggleAll( QUICKWIGGLE_ITERS )
        else
            structure.WiggleAll( QUICKWIGGLE_ITERS )
        end
        if math.abs( ss-getQuickScore( ) ) > scoreThreshold then 
            -- we've changed the structure "enough", send it off 
            break 
        end
    end
    SetCI( 1.0 )
    return math.abs( ss - getQuickScore( ) )
end

function SimpleWiggleShakeWrapper( minChangeWanted, pullingCI )
    if minChangeWanted == nil then minChangeWanted = BANDEDWIGGLE_TARGET end
    if pullingCI == nil then pullingCI = COREPULLING_CI end
    
    local scoreThreshold = math.min( minChangeWanted, getQuickScore( ) / 500.0)

    recentbest.Save( )
    local oldRbScore = getRBScore( )
    BandedWiggleUntilEnoughChange( scoreThreshold, pullingCI )
    local newRbScore = getRBScore( )
    if newRbScore > oldRbScore then recentbest.Restore( ) end

    DelBands( )
    ResetFrozenness( )
    recentbest.Save( )
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA ))
    recentbest.Restore( )
end

-- unused (BandedWiggleUntilEnoughChange is preferred)
function BandedWiggleUntilTargetMet( scoreDeltaTarget, maxTries, pullingCI )
    print( "BandedWiggleUntilTargetMet" )
    save.Quicksave( QS_ShortUseTemp1 )
    save.Quicksave( QS_ShortUseTemp2 )
    if pullingCI == nil then pullingCI = COREPULLING_CI end
    SetCI( pullingCI )

    local s1=getQuickScore( )
    local str=0.50
    local tries = 0
    local s2best = 99999
    while tries < maxTries do
        save.Quickload( QS_ShortUseTemp1)
        SetMyBandStrengths( str)
        structure.LocalWiggleAll( QUICKWIGGLE_ITERS ) 

        local s2 = math.abs( getQuickScore( ) - s1 )
        if s2 > scoreDeltaTarget * 1.2 then 
            str = str * ( 0.85 + random( 0.10 ) )
        elseif s2 < scoreDeltaTarget * 0.8 then 
            str = str + ( 1.0 - str ) * ( 0.05 + random( 0.1 ) )
        else 
            return math.abs( getQuickScore( ) - s1 )   -- we're good enough
        end
        if math.abs( s2 - scoreDeltaTarget ) < s2best then
            save.Quicksave( QS_ShortUseTemp2 )
            s2best = math.abs( s2 - scoreDeltaTarget )
        end

        tries = tries + 1
    end
    save.Quickload( QS_ShortUseTemp2 )  -- the best we could come up with...
    return math.abs( getQuickScore( ) - s1 )  -- returns the actual delta in case caller cares
end

function RegionWiggleBeforeAndAfter( startIdx, endIdx )
    recentbest.Save( )
    local oldRbScore = getRBScore( )

    SetCI( COREPULLING_CI )
    SimpleRegionWiggle( startIdx, endIdx, false, 2 + random( 4 ))
    local newRbScore = getRBScore( )
    if newRbScore > oldRbScore then recentbest.Restore( ) end

    DelBands( )  -- if we put in bands, they go now
    ResetFrozenness( )  -- if we froze anything, clean that up now

    SetCI( 1.0 )
    if randomBool( ) then
        qStabWiggle( randomBool( PROB_QSTAB_EXTRA ))
    else
        recentbest.Save( )
        ShakeOrMutate( SHAKESIDECHAIN_ITERS )
        recentbest.Restore( )
        SimpleRegionWiggle( startIdx, endIdx, true, random( 4, 8 ))
        recentbest.Restore( )
        LongWiggle( )  
    end
    recentbest.Restore( )
end

function WiggleShakeCleanQStab( )
    SimpleWiggleShake( )

    SetCI( 1.0 )
    DelBands( )  -- if we put in bands, they go now
    ResetFrozenness( )  -- if we froze anything, clean that up now

    recentbest.Save( )
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA ))
    recentbest.Restore( )
end

function PickBeforeAndAfterWiggler( doVariableStrength )
    if USE_LOWCI_WIGGLES then
        ShakeWiggleCleanShakeWiggle( )
    elseif doVariableStrength then
        SimpleWiggleShakeWrapper( BANDEDWIGGLE_TARGET, COREPULLING_CI ) 
    else
        WiggleShakeCleanQStab( )
    end
end

function nudge( )
    SetCI( 0.3 )
    structure.WiggleSelected( 1 )
    SetCI( 1 )
    structure.WiggleSelected( 8 )
end
function WiggleAndShakeForRebuilder( initScore, maxFall )
    print( "WiggleAndShakeForRebuilder" )
     -- initScore is the score before our rebuild happened
    recentbest.Save( )

    scoreA = getQuickScore( )
    if initScore - scoreA > maxFall then
        ShakeOrMutate( 2*SHAKESIDECHAIN_ITERS )
    elseif initScore - scoreA > maxFall/2 then
        ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    end

    scoreA = getQuickScore( )
    if randomBool( PROB_LOCALWIGGLE ) then
        structure.LocalWiggleAll( 15 )
    else
        structure.WiggleAll( 15 )
    end
    recentbest.Restore( )

    scoreB = getQuickScore( )
    if scoreB - scoreA < 10 and initScore - scoreB > 30 then
        -- The wiggle got stuck: try a nudge to see if that helps
        nudge( )
        recentbest.Restore( )
        scoreB = getQuickScore( )
    end

    if initScore - scoreB <= 100.0 then
        ShakeOrMutate( SHAKESIDECHAIN_ITERS )
        structure.WiggleSelected( 10 )
        recentbest.Restore( )
        nudge( )
        recentbest.Restore( )
    end   -- else not worth the effort
end

--------------------------------------------------
-- FUSES AND OTHER SLOW WIGGLERS
--------------------------------------------------
function BlueFuse( )
    recentbest.Save( )

    SetCI( .05 )
    structure.ShakeSidechainsAll( SHAKESIDECHAIN_ITERS )
    SetCI( 1 )
    structure.WiggleAll( 5 )
    recentbest.Restore( )

    SetCI( .07 )
    structure.ShakeSidechainsAll( SHAKESIDECHAIN_ITERS )
    SetCI( 1 )
    structure.WiggleAll( 5 )
    recentbest.Restore( )

    SetCI( .3 )
    structure.WiggleAll( 8 )
    SetCI( 1 )
    structure.WiggleAll( 15 )
    recentbest.Restore( )
end

function IterativeWiggleLocalByChunk(startSeg, endSeg, minGain, chunkSize, addFreeze)
    recentbest.Save( )
    currOffset = 0
    local lastScore = getQuickScore( )
    gain = minGain
    while gain >= minGain do
        if chunkSize == 1 then currOffset = 0
        else currOffset = random( chunkSize - 1 )
        end

        ctBlocks = math.floor( (endSeg - startSeg + 1) / chunkSize )
        for idx = 0,ctBlocks do
            local score = getQuickScore( )
            startIdx = startSeg + currOffset + chunkSize*idx
            WiggleRange( startIdx, math.min( startIdx + chunkSize, endSeg ), addFreeze, false )
            recentbest.Restore( )
            if addFreeze then ResetFrozenness( ) end
        end
        local currScore = getQuickScore( )
        gain = currScore - lastScore
        lastScore = currScore
    end
    recentbest.Restore( )
    if addFreeze then ResetFrozenness( ) end
end

function BasicCutAndWiggleFuse( baseIters, shift, offset )
    if shift == nil then shift = 3 end
    if offset == nil then offset = MinUnlockedSeg else offset = MinUnlockedSeg + (offset % shift) end

    for i = offset, MaxUnlockedSeg-1, shift do
        structure.InsertCut( i )
    end
    WiggleWiggleSimpleFuse( 0.30, 1.0, baseIters, 2*baseIters )
    ResetCuts( )

    WiggleWiggleSimpleFuse( 0.30, 1.0, baseIters, 4*baseIters )

    SetCI( 1.0 )
end

function IterativeCutAndWiggle( baseIters, minGain, shift )
    local score = getQuickScore( )
    local gain = minGain
    save.Quicksave( QS_ShortUseTemp1 )
    recentbest.Save( )
    repeat
        save.Quickload( QS_ShortUseTemp1 )
        BasicCutAndWiggleFuse( baseIters, shift, random( shift ) )
        local newscore = getQuickScore( )
        gain = 0.0
        if newscore >= score then 
            save.Quicksave( QS_ShortUseTemp1 )
            gain = newscore - score
            score = newscore
        end
    until gain < minGain
    save.Quickload( QS_ShortUseTemp1 )
    recentbest.Restore( )
end


-- BASIC FUSE WIGGLERS
function WiggleShakeSimpleFuse( ci )
    SetCI( ci )
    structure.WiggleAll( QUICKWIGGLE_ITERS )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
end
function WiggleWiggleSimpleFuse( ci1, ci2, shortiter, longiter )
    recentbest.Save( )
    local oldRbScore = getRBScore( )
    SetCI( ci1 )
    structure.WiggleAll( shortiter )
    local newRbScore = getRBScore( )
    if newRbScore > oldRbScore then recentbest.Restore( ) end

    SetCI( ci2 )
    recentbest.Save( )
    structure.WiggleAll( longiter )  -- "long" wiggle
    recentbest.Restore( )
end
function ShakeWiggleFuse( ci1, ci2 )  -- aka Fuze1
    SetCI( ci1 )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    SetCI( ci2 )
    LongWiggle( )  
end
function WiggleWiggleWiggleFuse( ci1, ci2 ) -- aka Fuze2
    SetCI( ci1 )
    LongWiggle( )  
    SetCI( 1.00 )
    LongWiggle( )  
    SetCI( ci2 )
    LongWiggle( )  
end
function WiggleShakeWiggleFuse( )  -- aka FuzeEnd
    SetCI( 1.00 )
    LongWiggle( )  
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    LongWiggle( )  
end

----------------------------------------------------------------------
--
---- BANDER FUNCTIONS
--
----------------------------------------------------------------------
function EnableAndDisableCysteineBands( )
	for i=1, CysteineBandPairCount do
		bIdx = InitialBandCount + COUNT_CYSTEINE_BANDS_USED*(i-1)
		if bIdx > band.GetCount() then
			print("bIdx too large - will not enable or disable")
		else
			pair = cysteinePairs[i]
			bondExists = DoesCysteineBondExist( pair[1], pair[2] )
			for j=1, COUNT_CYSTEINE_BANDS_USED do
				if bondExists then
					band.Disable( bIdx + j )
				else
					band.Enable( bIdx + j )
				end
			end
		end
	end
end

function DelBands( alsoCysteine )
	if alsoCysteine == nil then alsoCysteine = false end
	
    if alsoCysteine and not PRESERVE_USER_BANDS then
		print("Wiping all bands")
        band.DeleteAll( )
    else
        local allBandCount = band.GetCount( )
		firstBandToDel = 1 + InitialBandCount
		if not alsoCysteine then 
			firstBandToDel = firstBandToDel + 
				COUNT_CYSTEINE_BANDS_USED*CysteineBandPairCount
		end
        for i=allBandCount, firstBandToDel, -1 do 
			-- this *must* go in reverse order
			band.Delete( i ) 
		end
    end
end

function SetMyBandStrengths( strength )
	allBandCount = band.GetCount( )
	firstActionBand = InitialBandCount + 
			COUNT_CYSTEINE_BANDS_USED*CysteineBandPairCount + 1
    for i=firstActionBand, allBandCount do
        band.SetStrength( i, strength )
    end
end

function BandBetweenSegsWithParameters( seg1, seg2, strength, goalLength, atom1, atom2 )
    if not ( IsUnlockedSeg( seg1 ) or IsUnlockedSeg( seg2 ) ) then 
		return false 
	end
    if atom1 == nil then atom1 = CENTER_CARBON end
    if atom2 == nil then atom2 = CENTER_CARBON end

    local bIdx = band.AddBetweenSegments( seg1, seg2, atom1, atom2 )
    if bIdx ~= band.GetCount( ) then 
        DebugPrint( "failed to add band from "..seg1.." to "..seg2)
        return false 
    end
    if bIdx <= InitialBandCount then 
		-- SNH: don't accidentally change user-supplied bands
		return true 
	end

    if goalLength ~= nil then 
        band.SetGoalLength( bIdx, goalLength ) 
    end
    if strength ~= nil and strength > 0.0 then 
        band.SetStrength( bIdx, strength )
    end
    return true
end

function BandInSpaceWithParameters( seg, segFini, segInit, rho, theta, phi, strength, goalLength )
    if not IsUnlockedSeg( seg ) then return false end

    local bIdx = band.Add( seg, segFini, segInit, rho, theta, phi )
    if bIdx ~= band.GetCount( ) then return 0 end
    if bIdx <= InitialBandCount+COUNT_CYSTEINE_BANDS_USED*CysteineBandPairCount then 
		-- don't change user-supplied bands
		-- don't expect to hit our cysteine ones, but just to be safe...
		return 0 
	end 
    if goalLength ~= nil then 
        band.SetGoalLength( bIdx, goalLength ) 
    end
    if strength ~= nil and strength ~= 0.0 then 
        band.SetStrength( bIdx, strength ) 
    end
    return bIdx
end

-- fool with this later, does not work
function BandInSpaceFromAtomWithParameters( 
            seg, segFini, segInit, 
            rho, theta, phi, 
            strength, goalLength, atom)
    if not IsUnlockedSeg( seg ) then return false end
    local atomXAxis = 1  -- atoms 1,2,3 are on the backbone.
    local atomYAxis = 2
    if atom == 1 then atomXAxis, atomYAxis = 2,3 end
    if atom == 2 then atomXAxis, atomYAxis = 1,3 end
    local bIdx = band.Add( seg, segFini, segInit, rho, theta, phi, atom, atomXAxis, atomYAxis )
    if bIdx ~= band.GetCount( ) then return 0 end
    if bIdx <= InitialBandCount+COUNT_CYSTEINE_BANDS_USED*CysteineBandPairCount then 
		-- don't change user-supplied bands
		-- don't expect to hit our cysteine ones, but just to be safe...
		return 0 
	end 
    if goalLength ~= nil then 
        band.SetGoalLength( bIdx, goalLength ) 
    end
    if strength ~= nil and strength ~= 0.0 then 
        band.SetStrength( bIdx, strength ) 
    end
    return bIdx
end

function ZLB( strength, lower, upper, bandGivenRange )
    if lower == nil then lower = MaxUnlockedSeg + 1 end
    if upper == nil then upper = MinUnlockedSeg - 1 end
    if bandGivenRange == nil then bandGivenRange = false end     -- default is exclude a range

    if bandGivenRange then
        for i=lower, upper do 
            if i == 1 then 
                BandInSpaceWithParameters( i, i+1, i+2, 0.01, 0.0, 0.0, strength, 0.01)
            elseif i == segCt then 
                BandInSpaceWithParameters( i, i-1, i-2, 0.01, 0.0, 0.0, strength, 0.01)
            else
                BandInSpaceWithParameters( i, i+1, i-1, 0.01, 0.0, 0.0, strength, 0.01)
            end
        end
    else
       for i=MinUnlockedSeg+1, MaxUnlockedSeg-1 do 
            if i < lower or i > upper then
                BandInSpaceWithParameters( i, i+1, i-1, 0.01, 0.0, 0.0, strength, 0.01)
            end
        end
        if lower > 1 then
            BandInSpaceWithParameters( 1, 2, 3, 0.01, 0.0, 0.0, strength, 0.01)
        end
        if upper < segCt then
            BandInSpaceWithParameters( segCt, segCt-1, segCt-2, 0.01, 0.0, 0.0, strength, 0.01)
        end
    end
end

----------------------------------------------------------------------
-- SMALL-SCALE BANDERS
-- try making a cysteine joiner
----------------------------------------------------------------------
function PutBandInSpace( idx, maxRho, strength )
    local idx2, idx3
    if idx < segCt then
        idx2 = idx + 1
    else
        idx2 = idx - 2
    end
    if idx > 1 then
        idx3 = idx - 1
    else
        idx3 = idx + 2
    end

    -- random point in sphere of radius maxRho
    local theta, phi = randomThetaPhi( )
    local rho = ( maxRho * random( )^(1/3) ) + 0.001
    return BandInSpaceWithParameters( idx, idx2, idx3, rho, theta, phi, strength )
end

function PutBandToRandomSeg( idx, strength, atom )
    needNew = true
    local failedTries = 0
    while ( needNew and failedTries < 20 ) do
        idx2 = randomSeg( )  -- ok if this one isn't movable (we assume idx is movable)
        if idx2 > idx + 1 or idx2 < idx - 1 then 
            needNew = false 
            break
        else
            failedTries = failedTries + 1
        end
    end
    if not needNew then
        return BandBetweenSegsWithParameters( idx, idx2, strength, nil, atom, nil )
    end
    return false
end

function PutRandomCompBands()
    if not BAND_ADD_RANDOMS then return end
    if randomBool( ) then
        local ctRandBands = random( RANDOMBAND_COUNT - 2, RANDOMBAND_COUNT + 2 )
        local strength = random(0.5, 1.5) * BAND_STRENGTH_DEFAULT
        PutGentleCompBands( ctRandBands, 14.0, 30.0, 5, strength, 0.7 )  -- 0.7 is compression
    end
end

----------------------------------------------------------------------
-- MID-SCALE BANDERS
----------------------------------------------------------------------
function PutGlycineHingeBands( glycineIdx, strength )
    if glycineIdx == 1 or glycineIdx == segCt then return false end
    changeSucceeded = false
    selection.DeselectAll ( )

   -- construct the rods for the hinge (freeze before and after the glycine)
    local sR, eR = GetRegionStartAndEnd( glycineIdx + 1 )
    if eR < segCt and eR - glycineIdx < 3 then eR = eR + 1 end
    selection.SelectRange ( glycineIdx + 1 , eR )

    local sL, eL = GetRegionStartAndEnd( glycineIdx - 1 )
    if eL > 1 and glycineIdx - eL < 3 then eL = eL - 1 end
    selection.SelectRange ( sL, glycineIdx - 1 )
    ChangedFrozenness = true
    freeze.FreezeSelected ( true, true )
    selection.DeselectAll ( )

    -- put bands connecting the two rods 
    local maxLen = math.min( eR - sR + 1, eL - sL + 1 )
    if maxLen > 4 or (glycineIdx >= 4 and glycineIdx <= segCt - 4) then 
        maxLen = 4   -- try very hard to have 4 instead of only 2 or 3, even if it means bands beyond end of rods
    elseif maxLen < 1 then 
        return changeFailed  -- should only happen if glycine at an end
    end
    for i=1, maxLen do
        if BandBetweenSegsWithParameters( 
                              math.max(1, glycineIdx - i), 
                              math.min(segCt, glycineIdx + i), 
                              strength ) 
        then 
            changeSucceeded = true 
        end
    end
    return changeSucceeded
end

function PutGeneralHingeBands( startHinge, endHinge, strength )
    if startHinge > endHinge then startHinge, endHinge = endHinge, startHinge end
    if startHinge == 1 or endHinge == segCt then return false end

    changeSucceeded = false
    selection.DeselectAll( )

    -- construct the rods for the hinge (freeze before and after the region)
    local sR, eR = GetRegionStartAndEnd( endHinge + 1 )
    if eR < segCt and eR - endHinge < 3 then eR = eR + 1 end
    selection.SelectRange( endHinge + 1 , eR )
    local sL, eL = GetRegionStartAndEnd( startHinge - 1 )
    if eL > 1 and startHinge - eL < 3 then eL = eL - 1 end
    selection.SelectRange( sL, startHinge - 1 )

    ChangedFrozenness = true
    freeze.FreezeSelected( true, true )
    selection.DeselectAll( )

    -- put bands connecting the two rods 
    local maxLen = math.min( eR - sR + 1, eL - sL + 1 )
    if maxLen > 4 then maxLen = 4      -- we will surely have it be at least 2 by choice of eR and eL
    elseif maxLen < 1 then 
        -- should only happen if startHinge or endHinge is too close to an end
        return changeFailed 
    end
    for i=1, maxLen do
        if BandBetweenSegsWithParameters( 
                             math.max(1, startHinge - i), 
                             math.min(segCt, endHinge + i), 
                             strength ) 
        then 
            changeSucceeded = true 
        end
    end
    if not changeSucceeded then
        ResetFrozenness( )
    end
    return changeSucceeded
end

function PutCysteineBands( idx1, idx2 )
	if idx1 < 1 or idx1 > segCt or 
	   idx2 < 1 or idx2 > segCt or 
	   idx1 == idx2 then
		return false
	end
    if structure.GetAminoAcid( idx1 ) ~= 'c' or 
       structure.GetAminoAcid( idx2 ) ~= 'c' 
    then
        return false
    end

	at1 = GetTipAtomOfSeg( idx1 )
	at2 = GetTipAtomOfSeg( idx2 )
	-- fancy is to put more than one band, going after the sulfur atom:
	-- Put bands here:
	-- atS1 -> atS2     dist=d1=DISTANCE_CYSTEINE_BOND
	-- atS1 -> atS2-1   dist=d2=math.sqrt(d1*d1+1.8*1.8) 
	-- atS1-1 -> atS2   dist=d2
	-- atS1-1 -> atS2-1 dist=d3=math.sqrt(d2*d2+1.8*1.8)

	atS1 = GetCysteineSulfurAtom( idx1 )
	atS2 = GetCysteineSulfurAtom( idx2 )
    strength = CYSTEINE_BAND_STRENGTH
	d1 = DISTANCE_CYSTEINE_BOND
	if COUNT_CYSTEINE_BANDS_USED == 1 then
		-- my old way of doing this
		return BandBetweenSegsWithParameters( idx1, idx2, strength, d1, at1, at2 )
	else
		d2 = math.sqrt(d1*d1 + SCALE_CYSTEINE_OFF_1*SCALE_CYSTEINE_OFF_1) 
		d3 = math.sqrt(d2*d2 + SCALE_CYSTEINE_OFF_1*SCALE_CYSTEINE_OFF_1) 
		BandBetweenSegsWithParameters( idx1, idx2, strength, d1, atS1, atS2 )
		BandBetweenSegsWithParameters( idx1, idx2, strength, d2, atS1, atS2-1 )
		if COUNT_CYSTEINE_BANDS_USED >= 3 then
			BandBetweenSegsWithParameters( idx1, idx2, strength, d2, atS1-1, atS2 )
		end
		if COUNT_CYSTEINE_BANDS_USED >= 4 then
			BandBetweenSegsWithParameters( idx1, idx2, strength, d3, atS1-1, atS2-1 )
		end
	end
	return true  -- just assume it worked!
end

function PutCysteinePairBands( )
	changeSucceeded = false
	CysteineBandPairCount = 0
	for i=1, #cysteinePairs do
		pair = cysteinePairs[i]
		chSuc = PutCysteineBands( pair[1], pair[2] )
		if chSuc then
			-- will fail if user already has one of these in place
			CysteineBandPairCount = CysteineBandPairCount + 1
		end
		changeSucceeded = changeSucceeded or chSuc
	end
	return changeSucceeded
end

function PutPushPullBands( badSeg, bandsWanted, center, strength, wantPushNotPull )
    changeSucceeded = false

    local okSegList = {}
    local SC = structure.GetDistance( badSeg, center )
    local SR

    for i=1, segCt do
        local RC=structure.GetDistance( i, center )
        SR=structure.GetDistance( i, badSeg )
        if SR<15
            and RC<SC-3
            and math.abs( i - badSeg ) > 5
            and math.abs( i - center ) > 3
        then
            okSegList[#okSegList + 1] = i
         end
    end
    if #okSegList == 0 then
        print( "No suitable segs for idx="..badSeg.." with center="..center)
        return false
    end
    randomizeIndexList( okSegList )

    local bandsAdded = 0
    for i=1, #okSegList do
        SR=structure.GetDistance( i, badSeg )
        if wantPushNotPull then
            SR=SR + 4
        else
            SR=SR - 4
        end
        if SR>20 then SR=20 end
        if SR<3 then SR=3 end

        if BandBetweenSegsWithParameters( badSeg, okSegList[i], strength, SR ) then 
            changeSucceeded = true 
            bandsAdded = bandsAdded + 1
        end
        if bandsAdded >= bandsWanted then break end
    end --  FOR(i)
    return changeSucceeded
end

function PutVoidCrusherBands( idx, strength ) --make band if found void in area of segment "idx"
    changeSucceeded = false
    local t={} -- store possible segments

    for b=1,segCt do --test all segments
        if math.abs(idx-b) >= 15 then    -- else idx,b too near each other count-wise
            local ab = structure.GetDistance(idx, b)
            if ab > 10.0 then --no void if less
                local void=true

                for c=1,segCt do --searching for any segment between them
                    if c ~= idx and c ~= b then   -- don't want c = idx or b
                        local ac = structure.GetDistance(idx,c)
                        local bc = structure.GetDistance(b,c)
                        if 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
                        end
                    end
                end  -- END for c

                if void==true then 
                     t[#t+1]={idx,b}
                end
            end  -- END if idx and b are spatially at least 10 apart
        end -- END if idx and b are count-wise at least 15 apart
    end -- LOOP over candidate b's

    if #t>0 then
        for i=1,#t do
            local s=t[i][1]
            local e=t[i][2]
            local d=structure.GetDistance( s, e )
            d = d - 3.0  -- try to compress this much...
            if d<3 then d=3 end
            if d>20 then d=20 end
            if BandBetweenSegsWithParameters( s, e, strength, d ) then changeSucceeded = true end
        end
    end

    return changeSucceeded
end

function PutHelixUniformDistanceBands( startSeg, endSeg, shift, strength )
    changeSucceeded = false
    sumDist = 0
    for i=startSeg, endSeg-shift do
        sumDist = sumDist + structure.GetDistance( i, i+shift )
    end
    length = sumDist/( endSeg - shift - startSeg + 1 )

    for i=startSeg,endSeg-shift do
        if BandBetweenSegsWithParameters( i, i+shift, strength, length ) then changeSucceeded = true end
    end
     return changeSucceeded
end

function PutHelixFixerBands( startSeg, endSeg, strength, do4 )
    changeSucceeded = false

    if do4 then
        for i=startSeg,endSeg-4 do
            if BandBetweenSegsWithParameters( i, i+4, strength, HELIX_DISTANCE_SKIP4 ) then changeSucceeded = true end
        end
    else
        for i=startSeg,endSeg-5 do
            if BandBetweenSegsWithParameters( i, i+5, strength, HELIX_DISTANCE_SKIP5 ) then changeSucceeded = true end
        end
    end
    return changeSucceeded
end

function PutHelixCompressorBands( startSeg, endSeg, shift, strength, stretch )
    changeSucceeded = false
    for i=startSeg,endSeg-shift do
        dist = structure.GetDistance( i, i+shift )
        if BandBetweenSegsWithParameters( i, i+shift, strength, dist*stretch ) then changeSucceeded = true end
     end
     return changeSucceeded
end

function PutRotatorBands( startSeg, endSeg, length, strength, wantClockwise )
    changeSucceeded = false
    for i = math.max( startSeg , 2 ), math.min( endSeg , segCt - 1 ) do
        if wantClockwise then
            if BandInSpaceWithParameters( i, i + 1, i - 1, length, math.pi/2, math.rad ( 315 ), strength ) then
                changeSucceeded = true
            end
        else
            if BandInSpaceWithParameters( i, i - 1, i + 1, length, math.pi/2, math.rad ( 315 ), strength ) then
                changeSucceeded = true
            end
        end
    end
    return changeSucceeded
end

function PutPusherBands( startSeg, endSeg, length, strength, doForward )
    changeSucceeded = false
    init = math.max( 1, startSeg - 1 )
    fini = math.min( segCt, endSeg + 1 )
    if doForward then 
        for i = init + 1, fini - 1 do
            if BandInSpaceWithParameters( i, fini, init, length, math.pi/2, 0.0, strength ) then
                changeSucceeded = true
            end
        end
    else
        for i = fini - 1, init + 1, -1 do
            if BandInSpaceWithParameters( i, fini, init, length, math.pi/2, 0.0, strength ) then
                changeSucceeded = true
            end
        end
    end
    return changeSucceeded
end

function PutMoverBands( startSeg, endSeg, length, strength )
    changeSucceeded = false
    if startSeg > 1 then 
        init = startSeg - 1
    else
        init = startSeg
    end
    if endSeg < segCt then
       fini = endSeg + 1
    else
       fini = endSeg
    end

    theta, phi = randomThetaPhi( )
    for i = init + 1, fini - 1 do
        if BandInSpaceWithParameters( i, fini, init, length, theta, phi, strength ) then
            changeSucceeded = true
        end
    end
    return changeSucceeded
end

-- TODO consider putting more bands in "between" segments, with lengths that increase
-- as they head towards the midpoint
function PutFlexerBands( startSeg, endSeg, length, strength )
    if endSeg - startSeg < 3 then return false end
    changeSucceeded = false
    if startSeg > 1 then
        init = startSeg - 1
        first = startSeg
    else
        init = startSeg
        first = startSeg + 1
    end
    if endSeg < segCt then
       fini = endSeg + 1
       last = endSeg
    else
       fini = endSeg
       last = endSeg - 1
    end
    midpt = (startSeg + endSeg + 1) / 2

    -- initial band points inwards
    if BandInSpaceWithParameters( first, fini, init, 1.0, math.pi/2.0, 0.0, strength ) then
        changeSucceeded = true
    end

    --final band points inwards, too
    if BandInSpaceWithParameters( last, init, fini, 1.0, math.pi/2.0, 0.0, strength ) then
        changeSucceeded = true
    end

    -- middle one is longer and points perpendicularly to init-fini line. Consider
    -- replacing with ramping-up lengths over all "interior" segments.
    theta, phi = randomThetaPhi( )    -- ignoring phi, because we just want to use 0
    if BandInSpaceWithParameters( midpt, fini, init, length, theta, 0.0, strength ) then
        changeSucceeded = true
    end

    return changeSucceeded
end

function PutSheetStraightenerBands( startSeg, endSeg, strength )
    changeSucceeded = false
    for i=startSeg, endSeg-2 do
        changeSucceeded = BandBetweenSegsWithParameters( i, i+2, 
                                DISTANCE_INSHEET_N_TO_NPLUS2, strength ) 
                or changeSucceeded
    end
    return changeSucceeded
end

-- bandToNearest=false means "band to farthest"
function PutManyBandsToSeg( idx, maxBandCount, strength, compression, bandToNearest )
    changeSucceeded = false
    local segsToBand = {}
    local idxList = {}
    local distList = {}
    local totalDistance = 0.0
    local dist = 0.0

    -- we will construct segsToBand first, summing totalDistance as we go
    for i = 3, segCt - 3 do
        if i ~= idx and IsUnlockedSeg( i ) then
            if bandToNearest then
                if IsLocalMinimum( idx, i ) then
                    segsToBand[ #segsToBand + 1 ] = i
                    totalDistance = totalDistance + structure.GetDistance( idx, i )
                end
            else
                if IsLocalMaximum( idx, i ) then
                    segsToBand[ #segsToBand + 1 ] = i
                    totalDistance = totalDistance + structure.GetDistance( idx, i )
                end
            end
        end
    end

    -- now sort them so we can take only the "best" segments...
    for i=1, #segsToBand do idxList[ i ] = i end
    for i=1, #segsToBand do distList[ i ] = structure.GetDistance( idx, i ) end
    sortAllByScore( idxList, distList )

    -- now put the bands, preferring the "nearest" ones
    for i=1, math.min( #segsToBand, maxBandCount ) do
        local dist = structure.GetDistance( idx, segsToBand[ idxList[ i ] ] )
        if BandBetweenSegsWithParameters( idx, segsToBand[ idxList[ i ] ], strength, compression * dist ) then
            changeSucceeded = true
        end
    end

    return changeSucceeded
end


function PutTailDetacherBands( fromStart, length )
    if fromStart then
        s = 1
        e = 1 + length
        rel = 2 + length
    else
        s = segCt - length - 1
        e = segCt
        rel = segCt - length - 2
    end
    if rel < 1 or rel > segCt then return 0 end
    goalLen = 0
    phi = random( 0.5 * math.pi )
    theta = math.pi/2.0 + random( -0.5 * math.pi, 0.5 * math.pi )
    if fromStart then
        for i=e,s+1,-1 do
            goalLen = goalLen + 3.5
            BandInSpaceWithParameters( i, e+1, e+2, length, theta, phi, 1.0, goalLen ) 
         end
    else
        goalLen = goalLen + 3.5
        for i=s,e do
            BandInSpaceWithParameters( i, s-1, s-2, length, theta, phi, 1.0, goalLen )
        end
    end
    return true
end
function PutTailAttacherBands( fromStart, length )
    if fromStart then
        s = 1
        e = 1 + length
    else
        s = segCt - length - 1
        e = segCt
    end

    -- Idea is to find hydrophobics between s and e and attach bands
    -- from some of them to a "randomly" chosen nearby hydrophobic 
    -- not between s and e.
    local hookAas = { }
    for i=s, e do
        if structure.IsHydrophobic( i ) then hookAas[ #hookAas + 1] = i end
    end
    if  #hookAas == 0 then return false end  -- if we have no hydrophobics then give up for now
    randomizeIndexList( hookAas)
    local ctToHook = random(1, #hookAas )
    for i=1, ctToHook do
        local nearPhobics = {}
        toHook = hookAas[ i ]
        GetNearestPhobicsOutsideRange( toHook, s, e, nearPhobics, 5 )
        which = random( 1, #nearPhobics )
        BandBetweenSegsWithParameters( toHook, nearPhobics[which], 1.0, 2.0, BETA_CARBON, BETA_CARBON ) 
    end
    return true
end

----------------------------------------------------------------------
-- LARGE-SCALE BANDERS
----------------------------------------------------------------------
function PutSheetSewerBandsAllOver( maxDist, strength )
    changeSucceeded = false
    for seg1 = 2, segCt - 4 do
        for seg2 = seg1 + 3, segCt do
            if IsUnlockedSeg( seg1 ) or IsUnlockedSeg( seg2 ) then
                if structure.GetDistance( seg1, seg2 ) < maxDist and 
                   structure.GetSecondaryStructure( seg1 ) == 'E' and
                   structure.GetSecondaryStructure( seg2 ) == 'E'
                then
                    local changeSuc = BandBetweenSegsWithParameters(seg1, seg2, strength )
                    changeSucceeded = changeSucceeded or changeSuc
                end
            end
        end
    end
    return changeSucceeded
end

function PutCoreCompressBands(coreList, skipCt, startOffset, strength, compression, pullPhilics, doSidechains )
    if doSidechains == nil then doSidechains = false end
    if pullPhilics == nil then pullPhilics = false end

    local maxToBand = #coreList
    local bandList = {}
    local added = 0

    -- first put bands between all of our central segs
    -- notice: there's no protection from some of them being adjacent, but it's ok
    -- if one or more of these bands fail to be created.
    for i=1, maxToBand-1 do
        for j=i+1, maxToBand do
            if coreList[i] < coreList[j] - 1 or coreList[i] > coreList[j] + 1 then
                local goalLength= compression * structure.GetDistance( coreList[ i ], coreList[ j ] )
                if BandBetweenSegsWithParameters( coreList[ i ], coreList[ j ], strength, goalLength ) then
                    added = added + 1
                end
            end
        end
    end

    -- now put bands between the "rest" of the protein and our chosen ones
    for i=startOffset, segCt, skipCt do
        local skipThisI = false
        for j=1, maxToBand do 
            if coreList [ j ] == i then 
                skipThisI = true
                break
            end
        end

        if not skipThisI then 
            for j=1, maxToBand do
                if i < coreList[ j ] - 1 or i > coreList[ j ] + 1 then
                    goalLength= compression * structure.GetDistance( i, coreList[ j ] )
                    if pullPhilics and not structure.IsHydrophobic( i ) then
                        goalLength= goalLength + 1.0
                    end
                    atom = nil
                    if BAND_TO_SIDECHAINS and doSidechains then 
                        atom = PickAtomNumberForBand( i )
                    end
                    if BandBetweenSegsWithParameters( i, coreList[ j ], strength, goalLength, atom, nil ) then
                        added = added + 1
                    end
                end
            end
        end
    end

    return added > 0
end

function PutBandsEveryXToEveryY( Xrange, Xinit, Yrange, Yinit, strength, compression )
    local length = 0
    local added = 0
    for i=Xinit, segCt, Xrange do
        for j=Yinit, segCt, Yrange do
            if (i < j - 1 or i > j + 1) then
                if compression == nil then length = nil
                else length = compression * structure.GetDistance( i, j ) 
                end

                if BandBetweenSegsWithParameters(  i, j, strength, length ) then
                    added = added + 1
                end
            end
        end
    end
    return added > 0
end

function PutGentleCompBands( bandsWanted, minLength, maxLength, minSeparation, strength, compression )
    local bandsAdded = 0
    local failCount = 0
    repeat
        local s = randomMovableSeg( )
        local e = randomMovableFarSeg( s, minSeparation, minLength, maxLength )
        if e == 0 then 
            failCount = failCount + 1
        else
            -- yes, this would try to add the same segment twice, but that will fail simply...
            local dist = structure.GetDistance( s, e )
            changeSucceeded = BandBetweenSegsWithParameters( s,e, strength, dist * compression )
            if changeSucceeded then 
                bandsAdded = bandsAdded + 1
            else
                failCount = failCount + 1
            end
        end
    until bandsAdded == bandsWanted or failCount >= 100
    return true
end

----------------------------------------------------------------------
-- CONTACT-MAP BANDERS
----------------------------------------------------------------------
function PutAllContactMapBands( heatThreshold, madeStrength, madeLengthScale,
                                unmadeStrength, unmadeLengthScale, 
                                doSidechains, ditherValues, omitLongBands )
    changeSucceeded = false
    local baseStrength = strength
    for i = 1, segCt-2 do
        for j = i+2, segCt do
            if contactmap.GetHeat(i, j) >= heatThreshold then
                local goalLength = structure.GetDistance( i, j )
                if contactmap.IsContact( i, j ) then 
                    goalLength = madeLengthScale * goalLength
                    strength = madeStrength
                else
                    goalLength = unmadeLengthScale * goalLength
                    strength = unmadeStrength
                end
                if not omitLongBands or goalLength < MAX_ALLOWED_NONCONTACT then
                    if ditherValues then
                        goalLength = random( 0.95, 1.05, true ) * goalLength
                        strength = random( 0.95, 1.05, true ) * strength
                    end
                    local atom1 = nil
                    local atom2 = nil
                    if BAND_TO_SIDECHAINS and doSidechains then
                        atom1 = BETA_CARBON
                        atom2 = BETA_CARBON
                    end
                    changeSucceeded = BandBetweenSegsWithParameters( i, j, strength, goalLength, atom1, atom2 ) 
                                                  or changeSucceeded
                end
            end
        end
    end
    return changeSucceeded
end

function PutContactMapBands( )
    local threshold = CONTACTMAP_THRESHOLD
    local madeStrength = random( 0.96, 1.04, true ) * CONTACTMAP_CONTACTING_STRENGTH
    local madeLenScale = random( 0.96, 1.04, true ) * CONTACTMAP_CONTACTING_GOAL
    local unmadeStrength = random( 0.96, 1.04, true ) * CONTACTMAP_NONCONTACTING_STRENGTH
    local unmadeLenScale = random( 0.96, 1.04, true ) * CONTACTMAP_NONCONTACTING_GOAL
    local doSidechains =  BAND_TO_SIDECHAINS and randomBool( )
    local wantDither = randomBool( )
    local doVariableStrength = false
 
    return PutAllContactMapBands( threshold, madeStrength, madeLenScale,
                                                     unmadeStrength, unmadeLenScale, 
                                                     doSidechains, wantDither, true )
end

function PutSomeContactMapBands( heatThreshold, strength, ctBands, doSidechains, omitLongBands )
    changeSucceeded = false
    local hotList = {}
    for i = 1, segCt-2 do
        for j = i+2, segCt do
            local heat = contactmap.GetHeat(i, j)
            if heat >= heatThreshold and not contactmap.IsContact( i , j ) then
                if not omitLongBands or structure.GetDistance( i, j ) < MAX_ALLOWED_NONCONTACT then
                    hotList[ #hotList + 1] = { i, j, heat }
                end
            end
        end
    end

    randomizeIndexList( hotList )
    for i=1, math.min( ctBands, #hotList ) do
        local atom1 = nil
        local atom2 = nil
        if BAND_TO_SIDECHAINS and doSidechains then
            atom1 = PickAtomNumberForBand( hotList[i][1] )
            atom2 = PickAtomNumberForBand( hotList[i][2] )
        end
        local ch = BandBetweenSegsWithParameters( hotList[i][1], hotList[i][2], strength, nil, atom1, atom2 ) 
        changeSucceeded = ch or changeSucceeded
    end
    return changeSucceeded
end

function PutAllContactMapBandsToSeg( idx, heatThreshold, strength, doComp, doSidechain, omitLongBands )
    changeSucceeded = false
    local baseStrength = strength
    if doComp then
        -- COMP-style bands
        for i=1, segCt do
            if i < idx - 1 or i > idx + 1 then
                local heat = contactmap.GetHeat( i, idx )
                if heat >= heatThreshold then
                    local goalLength = structure.GetDistance( i, idx )
                    if contactmap.IsContact( i, idx ) then
                        strength = CONTACTMAP_CONTACTING_STRENGTH * baseStrength
                        goalLength = CONTACTMAP_CONTACTING_GOAL * goalLength
                    else
                        strength = CONTACTMAP_NONCONTACTING_STRENGTH * baseStrength
                        goalLength = CONTACTMAP_NONCONTACTING_GOAL * goalLength
                    end
                    if not omitLongBands or goalLength < MAX_ALLOWED_NONCONTACT then
                        local atom = nil
                        if BAND_TO_SIDECHAINS and doSidechain then
                            atom = PickAtomNumberForBand( i )
                        end
                        local changeSuc = BandBetweenSegsWithParameters( 
                                                      idx, i, strength, goalLength, 
                                                      nil, atom ) 
                        changeSucceeded = changeSuc or changeSucceeded
                    end
                end
            end
        end
    else
        -- only bands to non-contacting segs
        local hotList = {}
        for i = 1, segCt do
            if i < idx - 1 or i > idx + 1 then
                local heat = contactmap.GetHeat( i, idx )
                if heat >= heatThreshold and not contactmap.IsContact( i , idx ) then
                    if not omitLongBands or structure.GetDistance( i, idx ) < MAX_ALLOWED_NONCONTACT then
                        hotList[ #hotList + 1] = i
                    end
                end
            end
        end

        for i=1, #hotList do
            if hotList[ i ] ~= idx then  -- not supposed to happen!
                local atom = nil
                if BAND_TO_SIDECHAINS and doSidechain then
                    atom = PickAtomNumberForBand( hotList[i] )
                end
                local changeSuc = BandBetweenSegsWithParameters( 
                                                  idx, hotList[i], strength * contactmap.GetHeat( i, idx ), 
                                                  nil, nil, atom ) 
                changeSucceeded = changeSuc or changeSucceeded
            end
        end
    end
    return changeSucceeded
end

function PutContactMapBandsBetweenRegions( 
            heatThreshold, 
            startIdx1, endIdx1, startIdx2, endIdx2, 
            strength, doComp, doBetaCarbons, ditherValues, omitLongBands )
    changeSucceeded = false
    local baseStrength = strength
    local atom1 = CENTER_CARBON
    local atom2 = CENTER_CARBON
    if BAND_TO_SIDECHAINS and doBetaCarbons then
        atom1 = BETA_CARBON
        atom2 = BETA_CARBON
    end
    for i = startIdx1, endIdx1 do
        for j = startIdx2, endIdx2 do
            if contactmap.GetHeat( i, j ) >= heatThreshold then
                local goalLength = nil
                if doComp then
                    goalLength = structure.GetDistance( i, j )
                    if contactmap.IsContact( i, j ) then 
                        strength = CONTACTMAP_CONTACTING_STRENGTH * baseStrength
                        goalLength = CONTACTMAP_CONTACTING_GOAL * goalLength 
                    else
                        strength = CONTACTMAP_NONCONTACTING_STRENGTH * baseStrength
                        goalLength = CONTACTMAP_NONCONTACTING_GOAL* goalLength 
                    end
                end
                if not omitLongBands or goalLength < MAX_ALLOWED_NONCONTACT then
                    if ditherValues then
                        strength = random( 0.95, 1.05, true ) * strength
                        goalLength = random( 0.95, 1.05, true ) * goalLength 
                    end
                    changeSucceeded = BandBetweenSegsWithParameters( i, j, strength, goalLength, atom1, atom2 ) 
                                                  or changeSucceeded
                end
            end
        end
    end
    return changeSucceeded
end

function PutPushPullContactMapBands(
            badSeg, bandsWanted, 
            center, strength, heatThreshold, 
            wantPushNotPull )
    changeSucceeded = false

    local okSegList = {}
    local SC = structure.GetDistance( badSeg, center )
    local SR

    for i=1, segCt do
        local RC=structure.GetDistance( i, center )
        SR=structure.GetDistance( i, badSeg )
        if SR<15
            and RC<SC-3
            and math.abs( i - badSeg ) > 5
            and math.abs( i - center ) > 3
        then
            okSegList[#okSegList + 1] = i
        end
    end
    randomizeIndexList( okSegList )

    local bandsAdded = 0
    for i=1, #okSegList do
        SR=structure.GetDistance( i, badSeg )
        if wantPushNotPull then
            SR=SR + 4
        else
            SR=SR - 4
        end
        if SR>20 then SR=20 end
        if SR<3 then SR=3 end

        if BandBetweenSegsWithParameters( badSeg, okSegList[i], strength, SR ) then 
            changeSucceeded = true 
            bandsAdded = bandsAdded + 1
        end
        if bandsAdded >= bandsWanted then break end
    end --  FOR(i)
    return changeSucceeded
end

function PutVoidCrusherContactMapBands( idx, strength, heatThreshold )
    changeSucceeded = false
    local t={} -- store possible segments

    for b=1,segCt do --test all segments
        if math.abs(idx-b) >= 15 and contactmap.GetHeat( idx, b ) >= heatThreshold then
            local ab = structure.GetDistance(idx, b)
            if ab > 10.0 then --no void if less
                local void=true

                for c=1,segCt do --searching for any segment between them
                    if c ~= idx and c ~= b then   -- don't want c = idx or b
                        local ac = structure.GetDistance(idx,c)
                        local bc = structure.GetDistance(b,c)
                        if 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
                        end
                    end
                end  -- END for c

                if void==true then 
                     t[#t+1]={ idx, b }
                end
            end  -- END if idx and b are spatially at least 10 apart
        end -- END if idx and b are count-wise at least 15 apart
    end -- LOOP over candidate b's

    if #t>0 then
        for i=1,#t do
            local s=t[i][1]
            local e=t[i][2]
            local d=structure.GetDistance( s, e )
            d = d - 3.0  -- try to compress this much...
            if d<3 then d=3 end
            if d>20 then d=20 end
            if BandBetweenSegsWithParameters( s, e, strength, d ) then changeSucceeded = true end
        end
    end

    return changeSucceeded
end

--------------------------------------------------------------
--
---- MUTATOR FUNCTIONS
-- TODO: Consider adding a "very similar" mutator
--------------------------------------------------------------
-- this one doesn't care what amino used to be at idx
function MutateToAny( idx )
    local aa = structure.GetAminoAcid( idx )
    local aanew = aa
    local i = 0
    while ( aanew == aa ) do
        i = random( aaCount )
        aanew = AA( i )
    end
    print( "mutating " .. aa .. " at " .. idx .. " to " .. aanew )
    structure.SetAminoAcid( idx, aanew )
    return true
end

-- this one preserves hydrophob/hydrophil
function MutateToSimilar( idx )
    local aa = structure.GetAminoAcid( idx )
    local i = 0
    local aanew = aa
    while (aanew == aa) do
        if structure.IsHydrophobic( idx ) then
            i = random( aaPhobCount )
            aanew = AAphobic( i )
        else
            i = random( aaPhilCount )
            aanew = AAphilic( i )
        end
    end
    print( "mutating " .. aa .. " at " .. idx .. " to " .. aanew )
    structure.SetAminoAcid( idx, aanew )
    return true
end

----------------------------------------------------------------------
--
---- IDEALIZERS AND REBUILDERS
--
----------------------------------------------------------------------

----------------------------------------------------------------------
-- IDEALIZERS
----------------------------------------------------------------------

function IdealizeSelected( )
	recentbest.Save()
    band.DisableAll( )
    structure.IdealizeSelected( )
    band.EnableAll( )
	EnableAndDisableCysteineBands( )
	if not CheckCysteineSanity() then
		print("Rejected idealize due to cysteine insanity")
		recentbest.Restore()
		return false
	end
	return true
end

function IdealizeInRange( startSeg, endSeg, useCuts )
    if useCuts == nil then useCuts = true end
    local idxList = { }
    changeFailed = false
    if startSeg > endSeg then startSeg, endSeg = endSeg, startSeg end
    if startSeg < 1 then startSeg = 1 end
    if endSeg > segCt then endSeg = segCt end
    local tries = 0

    save.Quicksave( QS_ShortUseTemp1 )

    local delta = 0.0
    repeat
        save.Quickload( QS_ShortUseTemp1 )
        changeFailed = false
        if IDEALIZE_MAXTRIES > 1 and tries == IDEALIZE_MAXTRIES - 1 then
            delta = 0.0
            useCuts = false   -- overrule the user's request for use of cuts (should always succeed)
        end

        if useCuts then
            if startSeg - 2 > 1 then 
                structure.InsertCut( startSeg - 2 ) 
            end
            if endSeg + 1 < segCt then 
                structure.InsertCut( endSeg + 1 ) 
            end

            if startSeg - 3 >= 1 then
                BandBetweenSegsWithParameters( startSeg-3, startSeg, 1.0 + delta, 
                                               structure.GetDistance( startSeg-3, startSeg ) )
            end

            if endSeg + 3 <= segCt then
                BandBetweenSegsWithParameters( endSeg, endSeg + 3, 1.0 + delta, 
                                               structure.GetDistance( endSeg, endSeg+3 ) )
            end
        end

        selection.DeselectAll( )
        getSegsWithinRangeOfTube( idxList, startSeg, endSeg, IDEALIZE_TUBE_RADIUS )
        for i=1, #idxList do selection.Select( idxList[i] ) end
        changeSucceeded = IdealizeSelected( )
        selection.DeselectAll( )
		if changeSucceeded then
			if useCuts then
				-- get some wiggling happening
				WiggleZlbNicely( IDEALIZE_ZLB_STRENGTH + delta, 
										  math.max( 1, startSeg-2 ), math.min( segCt, endSeg+2 ) )

				-- now put the protein back in its proper shape
				if startSeg - 2 > 1 then
					if not DeleteCutAtIndex(startSeg - 2 ) and IDEALIZE_FORBID_UNHEALED_CUTS then 
						changeFailed = true 
					end
				end
				if changeFailed == false and endSeg + 1 < segCt then 
					if not DeleteCutAtIndex(endSeg + 1 ) and IDEALIZE_FORBID_UNHEALED_CUTS then 
						changeFailed = true 
					end
				end
			else    -- since we didn't put in any cuts, we can try simple wiggling to see how we did
				recentbest.Save( )
				SimpleRegionWiggle( startSeg, endSeg, randomBool() )
				ShakeOrMutate( SHAKESIDECHAIN_ITERS )
				recentbest.Restore( )
			end
		end

        if changeSucceeded and not changeFailed then 
            DelBands( )
            return true 
        end
        tries = tries + 1
        delta = delta + 0.25
    until tries >= IDEALIZE_MAXTRIES
    DelBands( )
    return false
end

function IdealizeInRange2( startSeg, endSeg )
    local idxList = { }
    save.Quicksave( QS_ShortUseTemp1 )

    if startSeg > 2 then structure.InsertCut( startSeg - 1 ) end
    if endSeg < segCt then structure.InsertCut( endSeg ) end

    selection.DeselectAll( )
    selection.SelectRange( startSeg, endSeg )       
    changeSucceeded = IdealizeSelected( )
	if changeSucceeded then
		if startSeg > 2 then 
			if not DeleteCutAtIndex( startSeg - 1 ) then
				changeFailed = true
			end
		end
		if endSeg < segCt then 
			if not DeleteCutAtIndex( endSeg ) then 
				changeFailed = true
			end
		end
	end
    if changeFailed or not changeSucceeded then
        save.Quickload( QS_ShortUseTemp1 )
        return false -- just give up
    end

    selection.DeselectAll( )
    recentbest.Save( )
    getSegsWithinRangeOfTube( idxList, startSeg, endSeg, IDEALIZE_TUBE_BIG_RADIUS )
    for i=1, #idxList do 
		selection.Select( idxList[i] ) 
	end
    structure.WiggleSelected( LONGWIGGLE_ITERS )

    selection.DeselectAll( )
    recentbest.Restore( )
	return true
end

function lowCiIdealizer()
    if IdealizeMin == nil then IdealizeMin = 1 end
    if IdealizeMax == nil then IdealizeMax = segCt end
    
    ExactSetCI( 0.05 )
    ZLB( IDEALIZEALL_ZLB_STRENGTH )
    selection.SelectRange( IdealizeMin, IdealizeMax )
    changeSucceeded = IdealizeSelected( )
	if changeSucceeded then
		structure.WiggleAll( 30 )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		structure.WiggleAll( 30 )
		selection.DeselectAll( )
	end
    DelBands( )
    print( "After Initial Idealizer = "..TrimNum( getScore( )))
	return changeSucceeded
end

----------------------------------------------------------------------
-- REBUILDERS
----------------------------------------------------------------------
function RebuildSelected( iters )
	recentbest.Save()
    band.DisableAll( )
    structure.RebuildSelected( iters )
    band.EnableAll( )
	EnableAndDisableCysteineBands( )
	if not CheckCysteineSanity() then
		print("Rejected rebuild due to cysteine insanity")
		recentbest.Restore()
		return false
	end
	return getQuickScore( ) ~= initscore 
end

function GetRebuildRange( idx, maxDist )
    if maxDist == nil then maxDist = math.floor( REBUILD_MAX_SEGLEN / 2 ) end
    
    local startIdx = math.max( 1, idx - random( maxDist ))
    local endIdx = math.min( segCt, idx + random( maxDist ))

    for i=idx, endIdx do
        if not IsMovableSeg(i) then 
            endIdx = i - 1
            break
        end
    end
    for i=idx, startIdx, -1 do
        if not IsMovableSeg(i) then 
            startIdx = i + 1
            break
        end
    end

    return startIdx, endIdx
end

function RebuildInRange( idxStart, idxEnd, maxTries, maxFall )
    if maxFall == nil then maxFall = REBUILD_MAX_MEDIUM_SCOREDROP end
    if maxTries == nil then maxTries = REBUILD_TRYCOUNT end
    wantHelix = isPossibleHelix( idxStart )
    
    changeSucceeded = false
    local iters = QUICKWIGGLE_ITERS
    local oldScoreMat = {}
    local locss = {}
    local changedSS = false
    local startScore = getQuickScore( )

    ConstructScoreMatrix( oldScoreMat )

    -- probably PickRebuild should decide whether to do secondary structure stuff
    for i=idxStart, idxEnd do locss[i] = structure.GetSecondaryStructure( i ) end
    for i=idxStart, idxEnd do
        if wantHelix then
            -- we think there's a helix, so tell rebuild to go with it
            changedSS = true 
            structure.SetSecondaryStructure( i, 'H' )
        elseif structure.GetSecondaryStructure( i ) == 'E' then 
            -- sheet isn't worth telling rebuild about
            changedSS = true 
            structure.SetSecondaryStructure( i, 'L' )
        end
    end

    selection.DeselectAll( )
    save.Quicksave( QS_ShortUseTemp1 )  -- holds our "running" pose
    save.Quicksave( QS_ShortUseTemp2 )  -- holds our "best so far" result

    bestRescore = -10000.0
    for i=1, maxTries do
        local idxListInner = { }
        save.Quickload( QS_ShortUseTemp1 )  -- bring up our "running" pose

        selection.SelectRange( idxStart, idxEnd )
        -- this param of "i" is weird: shouldn't this be iters?
        changeSucceeded = RebuildSelected( i )
        selection.DeselectAll( )

        if changeSucceeded then
              -- give the whole thing a little bit of shake-and-wiggle-time (deliberately not complete)
            structure.WiggleSelected( iters )
            getSegsWithinRangeOfTube( idxListInner, idxStart, idxEnd, REBUILD_INNER_TUBE_RADIUS )
            ShakeOrMutate( SHAKESIDECHAIN_ITERS, idxListInner )  -- only shake the "nearby" segments
            structure.WiggleAll( iters, true, true )  -- just a short bit of time on this wiggle to settle things down a bit

            rescore = getQuickScore( )
            if rescore > startScore then   -- we have a real improvement!
                save.Quicksave( QS_ShortUseTemp1 )  -- replace our "running" pose with this better one
                save.Quicksave( QS_ShortUseTemp2 ) -- our best-so-far is better
                startScore = rescore
                bestRescore = rescore
            end
        end
    end

    save.Quickload( QS_ShortUseTemp2 )
    if changedSS then
        for i=idxStart, idxEnd do structure.SetSecondaryStructure( i, locss[ i ] ) end
    end
    selection.DeselectAll( )
    return changeSucceeded
end

-- this does not have a quake, but it is the rebuild action seen in quake variants
function RebuildSection_FromQuaker( idxStart, idxEnd, maxTries, maxFall )
	local changeSucceeded = false
    local initScore = getQuickScore( )
    local tries = 0
    repeat
        selection.DeselectAll( )
        selection.SelectRange( idxStart, idxEnd )
        tries = tries + 1
        chSuc = RebuildSelected( tries ) -- again, steadily increasing iterations?
		changeSucceeded = changeSucceeded or chSuc
        newScore = getQuickScore( )
    until newScore ~= initScore or tries > maxTries

	if changeSucceeded then
		recentbest.Save( )   -- N.B. This is what recentbest was made for!
		RebuildSelected( maxTries - 1 )
		recentbest.Restore( )

		WiggleAndShakeForRebuilder( initScore, maxFall )
		-- quaking rebuild put a perform-quake-if-score-not-better-but-not-too-bad here
	end
    return changeSucceeded
end

function RebuildSimple( idxStart, idxEnd, addIdealize )
    if idxStart < 1 then idxStart = 1 end
    if idxEnd > segCt then idxEnd = segCt end
	if addIdealize == nil then addIdealize = REBUILD_ADDIDEALIZE end

    local initscore = getQuickScore( )
    local idxList = { idxStart }
    len = idxEnd - idxStart - 1
    for i=1, idxEnd - 1 do
        if idxStart + i <= segCt then idxList[ #idxList + 1 ] = idxStart + i
        else break
        end
    end

    selection.DeselectAll( )
    selection.SelectRange( idxStart, idxEnd )
    changeSucceeded = RebuildSelected( 1 )

	if changeSucceeded then
    print( "RS-Wiggler")
		recentbest.Save( )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS, idxList )
		recentbest.Restore( )
		if addIdealize then 
			chSuc = IdealizeSelected( )
			if chSuc then 
				recentbest.Restore( )
				ShakeOrMutate( SHAKESIDECHAIN_ITERS, idxList )
				recentbest.Restore( )
				structure.WiggleAll( 25 )
				recentbest.Restore( )
				ShakeOrMutate( SHAKESIDECHAIN_ITERS )
			end
			recentbest.Restore( )
		end
		LongWiggle( )
		recentbest.Restore( )
	end

    return changeSucceeded
end

function RebuildSimple2( idxStart, idxEnd )
    local initscore = getQuickScore( )
    selection.DeselectAll( )
    selection.SelectRange( idxStart, idxEnd)
    changeSucceeded = RebuildSelected( 1 )

	if changeSucceeded then
		print( "RS2-Wiggler")
		recentbest.Save( )
		SetCI( 0.5 )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		SetCI( 1.0 )
		structure.WiggleAll( 15 )
		recentbest.Restore( )
		SetCI( 0.5 )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		SetCI( 1.0 )
		LongWiggle( )
		recentbest.Restore( )
	end

    return changeSucceeded
end

function RebuildSimple3( idxStart, idxEnd, addIdealize )
	if addIdealize == nil then addIdealize = REBUILD_ADDIDEALIZE end
	
    local initscore = getQuickScore( )
	recentbest.Save()
    selection.DeselectAll( )
    selection.SelectRange( idxStart, idxEnd )
    changeSucceeded = RebuildSelected( 4 )

    if changeSucceeded then
		print( "RS3-Wiggler")
        recentbest.Save( )
        if addIdealize then 
            chSuc = IdealizeSelected( ) 
			if not chSuc then
				recentbest.Restore()
			end
        end
        VerySimpleShakeWiggle( )
        recentbest.Restore( )
    end
	return changeSucceeded
end

function RebuildBitspawn( idxStart, idxEnd, useLowCIRules )
    local initscore = getQuickScore( )
    selection.DeselectAll( )
    selection.SelectRange( idxStart, idxEnd )
    changeSucceeded = RebuildSelected( 4 )
    selection.DeselectAll( )

    if changeSucceeded then
		print( "RBits-Wiggler")
        if useLowCIRules then
            local ciOld = behavior.GetClashImportance( )
            ExactSetCI( 0.05 )
            ShakeOrMutate( SHAKESIDECHAIN_ITERS )
            ExactSetCI( ciOld )
            structure.WiggleAll( 16 )
        else
            ShakeOrMutate( SHAKESIDECHAIN_ITERS )
            structure.WiggleAll( 2 )
            ShakeOrMutate( SHAKESIDECHAIN_ITERS )
            structure.WiggleAll( 16 )
        end
        return true
    else
        return false
    end
	return changeSucceeded
end

function WiskyRebuildAllAction( )
    SelectAllMovable( )
    changeSucceeded = RebuildSelected( 1 )
    if changeSucceeded then
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		structure.WiggleAll( 10 )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		LongWiggle( )
	end
	return changeSucceeded
end

function RebuildAllWithIdealizeAction( )
    SelectAllMovable( )
    changeSucceeded = RebuildSelected( 1 )
    if changeSucceeded then
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		structure.WiggleAll( 10 )
		ShakeOrMutate( SHAKESIDECHAIN_ITERS )
		IterativeWiggleAll( true, true )
    
		lowCiIdealizer()
	end
	return changeSucceeded
end

----------------------------------------------------------------------
-- FUNCTIONS TO PERFORM SOME ACTUAL ACTION
----------------------------------------------------------------------
function PerformRebuildRangeAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildRange on range "..idxStart..","..idxEnd..
           " Idealize? "..BoolStr(REBUILD_ADDIDEALIZE) )

    local maxFall = REBUILD_MAX_MEDIUM_SCOREDROP
    changeSucceeded = false

    local initScore = getQuickScore( )
    changeSucceeded = RebuildInRange( idxStart, idxEnd, REBUILD_TRYCOUNT, maxFall )
    if changeSucceeded then
        if REBUILD_ADDIDEALIZE then
            save.Quicksave( QS_ShortUseTemp2 )  -- so we can throw out the idealize if we don't like result
            local score = getQuickScore( )
            local check = IdealizeInRange( idxStart, idxEnd, REBUILD_USECUTS )
            if check and score < getQuickScore( ) then 
                save.Quicksave( QS_ShortUseTemp2 )
                score = getQuickScore( )
            end
            save.Quickload( QS_ShortUseTemp2 )
        end

        local idxListOuter = {}
        getSegsWithinRangeOfTube( idxListOuter, idxStart, idxEnd, REBUILD_OUTER_TUBE_RADIUS )
        WiggleShakeWiggleList( idxListOuter  )
    end

    return changeSucceeded
end

function PerformRebuildSSAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    
    -- This tries to get a range matching some Secondary Structure:
    idxStart, idxEnd = GetRegionForOperation( idx )
    print( "RebuildSS on range " ..idxStart.. "," ..idxEnd..
           " Idealize? "..BoolStr(REBUILD_ADDIDEALIZE))

    changeSucceeded = false

    local initScore = getQuickScore( )
    changeSucceeded = RebuildInRange( idxStart, idxEnd, 
                                      REBUILD_TRYCOUNT, 
                                      REBUILD_MAX_MEDIUM_SCOREDROP )
    if changeSucceeded then
        if REBUILD_ADDIDEALIZE then
            -- use quickslot so we can throw out the idealize if we don't like result
            save.Quicksave( QS_ShortUseTemp2 )
            local preScore = getQuickScore( )
            local check = IdealizeInRange( idxStart, idxEnd, REBUILD_USECUTS )
            if (not check) or preScore < getQuickScore( ) then
                -- throw away the idealize result
                save.Quickload( QS_ShortUseTemp2 )
            end
        end

        local idxListOuter = {}
        getSegsWithinRangeOfTube( idxListOuter, idxStart, idxEnd,
                                  REBUILD_OUTER_TUBE_RADIUS )
        WiggleShakeWiggleList( idxListOuter )
    end

    return changeSucceeded
end

function PerformRebuildSimplePointAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    print( "RebuildSimple on seg "..idx )
    
    return RebuildSimple( idx, idx, REBUILD_ADDIDEALIZE )
end

function PerformRebuildSimpleRangeAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildSimpleRange on range " ..idxStart.. "," ..idxEnd )
    
    return RebuildSimple( idxStart, idxEnd, REBUILD_ADDIDEALIZE )
end

function PerformRebuildSimple2PointAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    print( "RebuildSimple2 on seg "..idx )
    
    return RebuildSimple2( idx, idx )
end

function PerformRebuildSimple2RangeAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildSimple2Range on range " ..idxStart.. "," ..idxEnd )
   
    return RebuildSimple2( idxStart, idxEnd )
end

function PerformRebuildSimple3RangeAction( idx ) 
    if idx == nil then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
     print( "RebuildSimple3Range on range " ..idxStart.. "," ..idxEnd )
   
    return RebuildSimple3( idxStart, idxEnd, randomBool( ))
end

function PerformRebuildFromQuakerAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildFromQuaker on range " ..idxStart.. "," ..idxEnd )
    
    return RebuildSection_FromQuaker( idxStart, idxEnd, 
                                      REBUILD_TRYCOUNT, 
                                      REBUILD_MAX_MEDIUM_SCOREDROP )
end

function PerformRebuildBitspawnAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildBitspawn on range " ..idxStart.. "," ..idxEnd )
   
    return RebuildBitspawn( idxStart, idxEnd, false )
end

function PerformRebuildBitspawnLowCI( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    idxStart, idxEnd = GetRebuildRange( idx )
    print( "RebuildBitspawnLowCI on range " ..idxStart.. "," ..idxEnd )
    
    return RebuildBitspawn( idxStart, idxEnd, true )
end

function PerformContactMapRebuildRangeAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local len = random( REBUILD_MAX_SEGLEN )
    if idx + len - 1 > segCt then len = segCt - idx + 1 end
    idxEnd = idx + len - 1
    print( "ContactMapRebuildRange on range ".. idx, idxEnd )

    changeSucceeded = PutContactMapBands( )
    if changeSucceeded then 
        if USE_LOWCI_WIGGLES then
            changeSucceeded = RebuildSimple3( idx, idxEnd, randomBool( ))
        elseif randomBool( ) then
            changeSucceeded = RebuildSimple( idx, idxEnd )
        else
            changeSucceeded = RebuildSimple2( idx, idxEnd )
        end
    end

    DelBands()
    qStabWiggle( randomBool(PROB_QSTAB_EXTRA)) 
    return changeSucceeded
end

-- think about all-loop and restoring secondary structure and such
function PerformContactMapRebuildSSAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local idxStart, idxEnd = GetRegionForOperation( idx )
    print( "ContactMapRebuildSS for range "..idxStart..","..idxEnd )

    local oldInnerTube = REBUILD_INNER_TUBE_RADIUS
    local oldOuterTube = REBUILD_OUTER_TUBE_RADIUS
    REBUILD_INNER_TUBE_RADIUS  = REBUILD_CM_INNER_TUBE_RADIUS
    REBUILD_OUTER_TUBE_RADIUS = REBUILD_CM_OUTER_TUBE_RADIUS

    changeSucceeded = PutContactMapBands( )
    if changeSucceeded then 
        if randomBool( ) then
            changeSucceeded = RebuildInRange( idxStart, idxEnd, 
                                              REBUILD_TRYCOUNT, 
                                              REBUILD_MAX_MEDIUM_SCOREDROP )
        else
            changeSucceeded = RebuildSection_FromQuaker( idxStart, idxEnd, 
                                                         REBUILD_TRYCOUNT, 
                                                         REBUILD_MAX_MEDIUM_SCOREDROP )
        end
    end

    REBUILD_INNER_TUBE_RADIUS = oldInnerTube
    REBUILD_OUTER_TUBE_RADIUS = oldOuterTube
    DelBands()
    qStabWiggle( randomBool(PROB_QSTAB_EXTRA)) 
    return changeSucceeded
end

function PerformContactMapRebuildAllAction( )
    print( "ContactMapRebuildAll" )

    changeSucceeded = PutContactMapBands( )
    if changeSucceeded then 
        WiskyRebuildAllAction( )
    end

    DelBands()
    qStabWiggle( true ) 
    return changeSucceeded
end

function PerformIdealizeRangeAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local maxSeg = math.min( segCt, idx + random( IDEALIZE_MAX_SEGLEN ) )
    local idxStart, idxEnd = GetRegionForOperation( idx )
    if idxEnd - idxStart + 1 > IDEALIZE_MAX_SEGLEN then 
         if idx - idxStart + 1 > IDEALIZE_MAX_SEGLEN then 
             idxStart = math.max( 1, idx - IDEALIZE_MAX_SEGLEN) 
        end
    end
    if idxEnd - idxStart + 1 > IDEALIZE_MAX_SEGLEN then 
         idxEnd = idxStart + IDEALIZE_MAX_SEGLEN - 1
    end
    print( "IdealizeRange for ".. idxStart..","..idxEnd )
    
    changeSucceeded = IdealizeInRange( idxStart, idxEnd, IDEALIZE_USE_CUTS )
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA) ) 

    return changeSucceeded
end

function PerformIdealizePointAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    print( "IdealizePoint for ".. idx )
    
    changeSucceeded = IdealizeInRange( idx, idx, IDEALIZE_USE_CUTS )
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA )) 

    return changeSucceeded
end

function PerformSingleBandInSpaceAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = random( BAND_STRENGTH_DEFAULT / 2, BAND_STRENGTH_DEFAULT, true )

    -- print is done here because PutBandInSpace is called by other pickers
    print( "BIS for segment " .. idx.." with strength "..TrimNum( strength ) )
    changeSucceeded = PutBandInSpace( idx, BIS_LENGTH, strength )
    if changeSucceeded then
        RegionWiggleBeforeAndAfter( idx, idx )
    end
    return changeSucceeded
end

function PerformBandToRandomSegAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_DEFAULT
    local atom = PickAtomNumberForBand( idx )
    local doVariableStrength = true

    -- print is done here because PutBandToRandomSeg is called by other pickers
    print( "BandBetween for segment " .. idx.." with strength "..TrimNum( strength ) )
    changeSucceeded = PutBandToRandomSeg( idx, strength, atom )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end


function PerformPushPullAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local segList = {}
    local wantPushNotPull = randomBool( )
    local ctBands = 5
    local strength = BAND_STRENGTH_DEFAULT
    local doVariableStrength = true

    print( "PushPull for segment " .. idx )
    getSegsOutsideSphere( segList, idx, 8 )
    if #segList == 0 then 
        -- nothing we can do!
        print("  PushPull failed due to lack of distant segs")
        return false 
    end
    local centerSeg = segList[ random( #segList) ]

    changeSucceeded = PutPushPullBands( idx, ctBands, centerSeg, strength, wantPushNotPull )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformVoidCrusherAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_DEFAULT
    local doVariableStrength = true

    print( "VoidCrusher for segment " .. idx )
    changeSucceeded = PutVoidCrusherBands( idx, strength )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformGlycineHingeAction( idx )
    if idx == nil then idx = randomAASeg( ) end
    if idx == 0 then return false end
    local strength = random( 0.25, BAND_STRENGTH_DEFAULT, true )

    print( "GlycineHinge for idx="..idx)
    
    local startSeg = math.max( 1, idx - 4 - random( 0, 2 ))
    local endSeg = math.min( segCt, idx + 4 + random( 0, 2 ))
    changeSucceeded = PutGlycineHingeBands( idx, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        -- had been SimpleRegionWiggle, but seemed to be doing nothing
        RegionWiggleBeforeAndAfter( startSeg, endSeg )
    end   
    return changeSucceeded
end

function PerformGeneralHingeAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local bisLen = BIS_LENGTH

    local startSeg, endSeg = GetRegionForOperation( idx )
    if endSeg - startSeg <= 1 then return false end -- something's very wrong
    if startSeg < idx - 2 then startSeg = idx - random( 0, 2 ) end
    if endSeg > idx + 2 then endSeg = idx + random( 0, 2 ) end
    print( "GeneralHinge for range "..startSeg..","..endSeg)

    changeSucceeded = PutGeneralHingeBands( startSeg, endSeg, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI ) 
        -- had been SimpleRegionWiggle, but seemed to be doing nothing
        RegionWiggleBeforeAndAfter( startSeg, endSeg )
    end    
    return changeSucceeded
end

function PerformRotatorAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local bisLen = BIS_LENGTH
    local goClockwise = randomBool()
    local doubleRotate = randomBool()

    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "Rotator for range "..startSeg..","..endSeg)
    if endSeg - startSeg <= 1 then return false end -- something's very wrong

    -- Put the bands, use wiggle to do a rotation
    changeSucceeded = PutRotatorBands( startSeg, endSeg, bisLen, strength, goClockwise )
    if changeSucceeded then
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
    end
    DelBands( )
    if doubleRotate then
        changeSucceeded = PutRotatorBands( startSeg, endSeg, bisLen, strength, goClockwise )
        if changeSucceeded then 
            SetCI( REGIONMOD_CI )
            SimpleRegionWiggle( startSeg, endSeg, false, 4 + random( 3 ) )
        end
    end
    DelBands( )
    CleanShakeQStab()
    return changeSucceeded
end

function PerformPusherAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local bisLen = BIS_LENGTH_LONG
    local goForward = randomBool()
    local doublePush = randomBool()

    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "Pusher for range "..startSeg..","..endSeg)
    if endSeg - startSeg <= 1 then return false end -- something's very wrong

    changeSucceeded = PutPusherBands( startSeg, endSeg, bisLen, strength, goForward )
    if changeSucceeded then
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
    end
    DelBands( )
    if doublePush then
        changeSucceeded = PutPusherBands( startSeg, endSeg, bisLen, strength, goForward )
        if changeSucceeded then 
            SetCI( REGIONMOD_CI )
            SimpleRegionWiggle( startSeg, endSeg, false, 4 + random( 3 ) )
        end
    end
    DelBands( )
    CleanShakeQStab()

    return changeSucceeded
end

function PerformMoverAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local bisLen = BIS_LENGTH_LONG
    local doubleMove = randomBool()

    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "Mover for range "..startSeg..","..endSeg)
    if endSeg - startSeg <= 1 then return false end -- something's very wrong

    changeSucceeded = PutMoverBands( startSeg, endSeg, bisLen, strength )
    if changeSucceeded then
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
    end
    DelBands( )
    if doubleMove then
        changeSucceeded = PutMoverBands( startSeg, endSeg, bisLen, strength )
        if changeSucceeded then 
            SetCI( REGIONMOD_CI )
            SimpleRegionWiggle( startSeg, endSeg, false, 4 + random( 3 ) )
        end
    end
    DelBands( )
    CleanShakeQStab()

    return changeSucceeded
end

function PerformFlexerAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local bisLen = BIS_LENGTH
    local doubleFlex = randomBool()

    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "Flexer for range "..startSeg..","..endSeg)
    if endSeg - startSeg <= 1 then return false end -- something's very wrong

    changeSucceeded = PutFlexerBands( startSeg, endSeg, bisLen, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
    end
    DelBands( )
    if doubleFlex then
        changeSucceeded = PutFlexerBands( startSeg, endSeg, bisLen, strength )
        if changeSucceeded then 
            SetCI( REGIONMOD_CI )
            SimpleRegionWiggle( startSeg, endSeg, false, 4 + random( 3 ) )
        end
    end
    DelBands( )
    CleanShakeQStab()

    return changeSucceeded
end

function PerformSheetStraightenerAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION

    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "SheetStraightener for range "..startSeg..","..endSeg)
    if endSeg - startSeg <= 1 then return false end -- something's very wrong

    changeSucceeded = PutSheetStraightenerBands( startSeg, endSeg, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
        DelBands( )
        CleanShakeQStab()
    end
    return changeSucceeded
end

function PerformHelixUniformAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "HelixUniform for range "..startSeg..","..endSeg)
    
    if endSeg - startSeg <= 1 then return false end -- something's very wrong
    local regionLength = endSeg - startSeg + 1
    local hShift = 4
    if randomBool( 0.25 ) then 
        hShift = math.min( regionLength - 2, randomDice( 4, 1, 2 ) ) -- a number 4 - 8, preferring 6ish
    end

    if startSeg > idx - 2 then startSeg = idx - random( 2, 4 ) end
    if endSeg < idx + 2 then endSeg = idx + random( 2, 4 ) end
    
    changeSucceeded = PutHelixUniformDistanceBands( startSeg, endSeg, hShift, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
        DelBands( )
        CleanShakeQStab()
    end
    return changeSucceeded
end

function PerformHelixCompressorAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local hComp = random( BAND_MINCOMPRESS, BAND_MAXEXPAND )
    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "HelixCompressor for range "..startSeg..","..endSeg)
    
    if endSeg - startSeg <= 1 then return false end -- something's very wrong
    local regionLength = endSeg - startSeg + 1
    local hShift = 4
    if randomBool( 0.25 ) then 
        hShift = math.min(regionLength - 2, randomDice(4, 1, 2)) -- a number 4 - 8, preferring 6ish
    end
    
    changeSucceeded = PutHelixCompressorBands( startSeg, endSeg, hShift, strength, hComp )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
        DelBands( )
        CleanShakeQStab()
    end
    return changeSucceeded
end

function PerformHelixFixerAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    local strength = BAND_STRENGTH_REGION
    local do4 = randomBool( )
    local startSeg, endSeg = GetRegionForOperation( idx )
    print( "HelixFixer for range "..startSeg..","..endSeg)
    
    if endSeg - startSeg <= 1 then return false end -- something's very wrong
    
    changeSucceeded = PutHelixFixerBands( startSeg, endSeg, strength, do4 )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        SimpleRegionWiggle( startSeg, endSeg, false, 2 + random( 4 ) )
        DelBands( )
        CleanShakeQStab()
    end
    return changeSucceeded
end

function PerformSheetSewerAction( )
    local maxSheetDistance = random( 4.0, 12.0, true )   -- default choice was 8.0
    local strength = BAND_STRENGTH_DEFAULT
    
    print( "SheetSewer")

    changeSucceeded = PutSheetSewerBandsAllOver( maxSheetDistance, strength )
    if changeSucceeded then 
        SetCI( REGIONMOD_CI )
        WiggleShakeCleanQStab( )
    end
    return changeSucceeded
end

function PerformMultiBanderAction( )
    local strength = random( 0, maxStrength, true ) + 0.001
    local ctBandsToMake = randomDice(3, 1, 4)  -- at least 3, maybe as many as 12
    local bandCt = 0
    local doVariableStrength = true
    
    print( "MultiBander with bands="..ctBandsToMake )
    while bandCt < ctBandsToMake do
        local which = random( 2 )
        if which == 1 then
            local seg = randomMovableSeg( )
            local maxRho = 1.0
            if PutBandInSpace( seg, maxRho, strength ) then
                changeSucceeded = true
                bandCt = bandCt + 1
            end
        else
            local idx = randomMovableSeg( )
            local atom = PickAtomNumberForBand( idx )
            if PutBandToRandomSeg( idx, strength, atom ) then
                bandCt = bandCt + 1
            end
        end
    end
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength )
    end
    return changeSucceeded
end

function PerformCentralBanderAction( )
    local strength = BAND_STRENGTH_DEFAULT
    local maxCenters = randomDice( 2, 1, 3)
    changeSucceeded = false
    local doVariableStrength = true

    local coreSegs = {}
    getCentralHydrophobics( coreSegs, maxCenters)
    local compression = randomDice( 2, BAND_MINCOMPRESS/2, BAND_MAXEXPAND/2 )
    local skipCt = random(4, 10)
    local startOffset = random(1, skipCt)

    print( "CentralBander with bands="..ctBandsToMake )
    changeSucceeded = PutCoreCompressBands( coreSegs, skipCt, startOffset, strength, 
                                               compression, randomBool( ), randomBool( ) )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformCoreStirrerAction( )
    -- not well suited to a lowCI situation
    changeSucceeded = false

    print( "CoreStirrer" )
    local coreSegs = {}
    getCentralHydrophobics( coreSegs, math.floor(segCt / 8) )   -- for a 100-seg protein, this is about 13

    -- perform an expansion
    local maxCenters = randomDice( 3, 1, 3 ) -- range 3 to 9, prefers 5-7
    local thisCore = {}
    local expansion = random( 1.00, BAND_MAXEXPAND )
    local strength = BAND_STRENGTH_DEFAULT
    local skipCt = random( 3, 13 )
    local startOffset = random( 1, skipCt )
    local pullPhilics = randomBool( )
    local doSidechains = randomBool( )

    getCentralHydrophobics( thisCore, maxCenters )
    changeSucceeded = PutCoreCompressBands( thisCore, skipCt, startOffset, strength, 
                                               expansion, pullPhilics, doSidechains )
    if not changeSucceeded then return false end

    print( "BandedWiggleUntilEnoughChange" )
    BandedWiggleUntilEnoughChange( getQuickScore( ) / 40.0, COREPULLING_CI )
    DelBands( )

    --
    -- jiggle a bunch of rotamers, chosen from phobics "near" the center.
    --
    local ctFlips = 0
    for i=1, #coreSegs do
        local ctRots = rotamer.GetCount( coreSegs[ i ] )
        if ctRots > 1 then 
            rotamer.SetRotamer( coreSegs[ i ], random( ctRots ) ) 
        end
    end

    qStabWiggle( true ) -- complete action of SimpleWiggleShakeWrapper now!
    return true
end

-- THIS one is normal band operation: deletes bands then wiggles out
function PerformEveryXToEveryYBandAction( )
    local compr
    compr = random( BAND_MINCOMPRESS, BAND_MAXEXPAND, true )
    local Xinterval = random( math.floor( segCt / 12), math.floor( segCt / 4 ) )
    local Yinterval = random( 2 + Xinterval, 2 + math.floor( segCt / 4 ) )
    local Xinit = random( 2*Xinterval )
    local Yinit = random( 2*Yinterval )
    local strength = random( 0.01, BAND_STRENGTH_DEFAULT, true )

    print( "BandEveryXToY for X,Y=" ..Xinterval.. "," ..Yinterval..
           " starting " ..Xinit.. "," ..Yinit )
    changeSucceeded = PutBandsEveryXToEveryY( Xinterval , Xinit, Yinterval, Yinit, strength, compr )
    if changeSucceeded then 
        BandedWiggleUntilEnoughChange( BANDEDWIGGLE_TARGET, COREPULLING_CI )
    end
    DelBands( )
    qStabWiggle( randomBool( PROB_QSTAB_EXTRA )) 
    return changeSucceeded
end

-- THIS ONE IS SPECIAL: it is called whenever we want to pre-expand our protein
function PerformEveryXToEveryYExpansion( )
    local compr = random( BAND_MINEXPAND, BAND_MAXEXPAND, true )
    local Xinterval = random( math.floor( segCt / 12), math.floor( segCt / 4 ) )
    local Yinterval = random( 2 + Xinterval, 2 + math.floor( segCt / 4 ) )
    local Xinit = random( 2*Xinterval )
    local Yinit = random( 2*Yinterval )
    local strength = random( 0.5, 1.50, true )  -- FIND OUT: is always strong good?
    print("XToYExpander: " ..Xinterval.. "x" ..Yinterval..
          " expand="..TrimNum(compr)..
          " strength="..TrimNum(strength))

    local initScore = getScore()
    changeSucceeded = PutBandsEveryXToEveryY( Xinterval , Xinit, Yinterval, Yinit, strength, compr )
    if changeSucceeded then 
        SimpleWiggleShake( )
    end
    DelBands( )
    local endScore = getScore()
    print("Expander Score drop = "..TrimNum(endScore-initScore))
    return changeSucceeded
end

function PerformGentleCompAction( )
    local minLength = random( 9.0, 16.0, true )      -- 10 and 14 seen in scripts
    local maxLength = random( 18.0, 32.0, true )   -- 20 and 30 seen in scripts; tempting to do ~2*minLength
    local minSeparation = random( 4, 12 )              -- 5 and 10 seen in scripts
    local compression = random( BAND_MINCOMPRESS, 0.90, true )  -- 0.70 seen in scripts
    local bandsWanted  = randomDice( 3, 2, 5 )     -- usually tunable; 7-9 for standalone is typical; 3-5 when added
    local strength = random( 0.50, 1.50) * BAND_STRENGTH_DEFAULT
    local doVariableStrength = true
    
    print( "GentleComp for " ..bandsWanted.. " bands" )
    changeSucceeded = PutGentleCompBands( bandsWanted, minLength, maxLength, 
                                          minSeparation, strength, compression )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformGentleCompLowCIAction( )
    local minLength = random( 9.0, 16.0, true )    -- 10 and 14 seen in scripts
    local maxLength = random( 18.0, 32.0, true )   -- 20 and 30 seen in scripts; tempting to do ~2*minLength
    local minSeparation = 15
    local compression = 1.0
    if randomBool() then 
        compression = random( BAND_MINCOMPRESS, 0.90, true )
    else
        compression = random( 1.10, BAND_MAXEXPAND, true )
    end
    local bandsWanted  = 3
    local strength = BAND_STRENGTH_DEFAULT
    local doVariableStrength = true
    
    print( "GentleCompLowCI for " ..bandsWanted.. " bands" )
    changeSucceeded = PutGentleCompBands( bandsWanted, minLength, maxLength, 
                                          minSeparation, strength, compression )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformLSQuakeAction( )
    local idx = randomSeg( )    -- we don't care if it's movable
    local wantFreeze = randomBool( )
    local ctBands = randomDice( 2, 2, 4 )   -- was 6
    local strength = BAND_STRENGTH_DEFAULT
    local compression = random( BAND_MINCOMPRESS, BAND_MAXEXPAND )  -- not present before. strength was tuned
    local doVariableStrength = true
    
    ResetFrozenness( )

    local startIdx, endIdx = GetRegionForOperation( idx )
    print( "LSQuaker over range "..startIdx..", "..endIdx )

    for i=startIdx, endIdx do
        local changeSuc =  PutManyBandsToSeg( i, ctBands, strength, compression, true)
        changeSucceeded = changeSucceeded or changeSuc
    end

    if changeSucceeded then
        if wantFreeze then 
            selection.DeselectAll( )
            selection.SelectRange( startIdx, endIdx )
            ChangedFrozenness = true
            freeze.FreezeSelected( true, true )
            selection.DeselectAll( )
        end
        PickBeforeAndAfterWiggler( doVariableStrength )
    end
    return changeSucceeded
end

function PerformTailWhipAction( idx )
    if idx == nil then fromStart = randomBool( ) end
    if idx < segCt / 2 then fromStart = true else fromStart = false end
    
    length = 0
    local r = randomBool()
    if r then
        if fromStart then idx = 1 else idx = segCt end
        st,ed = GetRegionForOperation( idx )
        length = ed - st + 1
    end
    if length < 5 then
        length = random( 5,15 )
    end

    if fromStart then
        str = "start"
    else
        str = "end"
    end
    print( "TailWhip at " .. str .. " length="..length )

    selection.DeselectAll( )
    
    s = 1
    e = length
    if fromStart then
        changeSucceeded = PutTailDetacherBands( true, length )
    else
        s = segCt - length + 1
        e = segCt
        changeSucceeded = PutTailDetacherBands( false, length )
    end
    selection.SelectRange( s, e )
    if not changeSucceeded then return false end
    SimpleRegionWiggle( s, e, false, random(4) )  -- pull the tail away from the main body
    DelBands( )
    selection.DeselectAll( )

    if fromStart then
        changeSucceeded = PutTailAttacherBands( true, length )
    else
        changeSucceeded = PutTailAttacherBands( false, length )
    end
    selection.SelectRange( s, e )
    if not changeSucceeded then return false end
    SimpleRegionWiggle( s, e, false, random(4) )  -- bring the tail back to the main body
    selection.DeselectAll( )
    DelBands( )

    -- do our own "after" wiggle
    recentbest.Save( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    IterativeWiggleAll( true, true )
    recentbest.Restore( )

    return true
end

function PerformSingleContactMapAction( idx )
    if idx == nil then idx = randomContactableSeg( ) end
    print( "SingleContactMapBand at ".. idx )
    changeSucceeded = false
    local heatThreshold = CONTACTMAP_THRESHOLD
    local strength = 1.5 * BAND_STRENGTH_DEFAULT
    local doVariableStrength = true

    changeSucceeded = PutSomeContactMapBands( heatThreshold, strength, 1, randomBool( ), true )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformSomeContactMapBandsAction( )
    print( "SomeContactMapBands" )
    changeSucceeded = false
    local heatThreshold = CONTACTMAP_THRESHOLD
    local strength = BAND_STRENGTH_DEFAULT
    local ctWanted = random( 5, 10 )
    local doVariableStrength = true

    changeSucceeded = PutSomeContactMapBands( heatThreshold, strength, ctWanted, randomBool( ), true )
    if changeSucceeded then 
        PutRandomCompBands() 
        PickBeforeAndAfterWiggler( doVariableStrength )
    end 
    return changeSucceeded
end

function PerformContactMapBandsForSegAction( idx )
    print( "ContactMapBandsForSeg for idx " ..idx )
    if idx == nil then idx = randomContactableSeg( ) end
    changeSucceeded = false
    local heatThreshold = CONTACTMAP_THRESHOLD
    local strength = BAND_STRENGTH_DEFAULT
    local doSidechain  = false
    local doVariableStrength = false
    if BAND_TO_SIDECHAINS then
        doSidechain = randomBool( )
    end

    changeSucceeded = PutAllContactMapBandsToSeg( idx, heatThreshold, strength, true, doSidechain, true )
    if changeSucceeded then 
        PutRandomCompBands() 
        PickBeforeAndAfterWiggler( doVariableStrength )
    end 
    return changeSucceeded
end

function PerformCMTwoRegionBandAction( idx1 )
    print( "CMTwoRegionBander for idx " ..idx1 )
    if idx1 == nil then idx1 = randomContactableSeg( ) end
    local strength = BAND_STRENGTH_DEFAULT
    local heatThreshold = CONTACTMAP_THRESHOLD
    local doVariableStrength = false

    local startIdx1, endIdx1 = GetRegionStartAndEnd( idx1 ) 
    if startIdx1 == 1 and endIdx1 == segCt then return false end -- we just can't deal with this
    local idx2 = startIdx1
    local tryCount = 0
    while idx2 >= startIdx1 and idx2 <= endIdx1 and tryCount < 50 do
        for i=startIdx1, endIdx1 do
            idx2 = randomContactSeg( i, CONTACTMAP_THRESHOLD )
            if idx2 ~= 0 and (idx2 < startIdx1 or idx2 > endIdx1) then break end -- got one
        end
        tryCount = tryCount + 1
    end
    if idx2 == 0 then return false end -- we have failed
    if tryCount == 50 then return false end
    local startIdx2, endIdx2 = GetRegionStartAndEnd( idx2 )

    changeSucceeded = PutContactMapBandsBetweenRegions( heatThreshold, 
                                                startIdx1, endIdx1, startIdx2, endIdx2, 
                                                strength, true, randomBool( ), randomBool( ), 
                                                true )
    if changeSucceeded then 
        PutRandomCompBands() 
        PickBeforeAndAfterWiggler( doVariableStrength )
    end 
    return changeSucceeded
end

function PerformContactMapperAction( )
    print("ContactMapper" )
    changeSucceeded = PutContactMapBands( )
    if changeSucceeded then 
        PutRandomCompBands() 
        PickBeforeAndAfterWiggler( doVariableStrength )
    end 
    return changeSucceeded
end

function PerformPushPullContactMapAction( badSeg )
    if badSeg == nil then badSeg = randomContactableSeg( ) end
    local segList = {}
    local wantPushNotPull = randomBool( )
    local ctBands = 5
    local strength = BAND_STRENGTH_DEFAULT
    local doVariableStrength = true

    print("PushPullCM at seg "..badSeg )
    getCMSegsOutsideSphere( segList, badSeg, 8, CONTACTMAP_THRESHOLD )
    if #segList == 0 then 
        -- nothing we can do!
        print( "  PushPull failed due to lack of distant 'center' segs" )
        return false 
    end
    centerSeg = segList[ random( #segList ) ]

    changeSucceeded = PutPushPullContactMapBands( badSeg, ctBands, centerSeg, strength, 
                                                             CONTACTMAP_THRESHOLD, wantPushNotPull )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

function PerformVoidCrusherContactMapAction( idx )
    if idx == nil then idx = randomContactableSeg( ) end
    local strength = BAND_STRENGTH_DEFAULT
    local doVariableStrength = true

    print("VoidCrusherCM at seg "..idx )
    changeSucceeded = PutVoidCrusherContactMapBands( idx, strength, CONTACTMAP_THRESHOLD )
    if changeSucceeded then 
        PickBeforeAndAfterWiggler( doVariableStrength ) 
    end
    return changeSucceeded
end

-- generally not used due to very expensive sequence of actions
function PerformLocalQuakeAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end
    print( "LocalQuake at ".. idx )
    local ctBands = randomDice( 2, 2, 4 )   -- was 6
    local compression = random( BAND_MINCOMPRESS, BAND_MAXEXPAND )
    local strength = BAND_STRENGTH_DEFAULT
    save.Quicksave( QS_ShortUseTemp1 )
    local bestScore = getQuickScore( )
    changeSucceeded = false
    local changeSuc = false

    local lqIdx = idx
    for i=1, 10 do
        save.Quickload( QS_ShortUseTemp1 )
        local currScore = getQuickScore( )
        DelBands( )
        changeSuc =  PutManyBandsToSeg( lqIdx, ctBands, strength, compression, true )
        if changeSuc then
            idx2 = math.min( segCt, lqIdx + 3 )
            changeSuc =  PutManyBandsToSeg( idx2, ctBands, strength, compression, true )
        end

        if changeSuc then
            if USE_LOWCI_WIGGLES then
                ShakeWiggleCleanShakeWiggle( )
            else
                SimpleWiggleShakeWrapper( BANDEDWIGGLE_TARGET, COREPULLING_CI )
            end
            newScore = getQuickScore( )
            if newScore > currScore then save.Quicksave( QS_ShortUseTemp1 ) end
            changeSucceeded = true  -- we only have to succeed once to count as a success
        end
        lqIdx = segCt - 2
        while lqIdx > segCt - 3 do
            lqIdx = randomSeg( )
        end
    end
    save.Quickload( QS_ShortUseTemp1 )
    return changeSucceeded
end

function PerformMutatorAction( idx )
    if idx == nil then idx = randomMovableSeg( ) end  -- better would be random mutable seg
    print( "SpotMutator for segment " .. idx )

    if randomBool( PROB_MUTATESIMILAR ) then
        changeSucceeded = MutateToSimilar( idx )
    else
        changeSucceeded = MutateToAny( idx )
    end
    PickBeforeAndAfterWiggler( false )
    return changeSucceeded
end


--------------------------------------------------------------------------------------------------------
-- NEIGHBOR-GENERATION
--------------------------------------------------------------------------------------------------------
function PickActionFromList( allowedActions, idx )
    local isRebuild = false
    local didExpand = false
	
    algName = ChooseAlgorithm( allowedActions )
	
 	if PERFORM_CYSTEINE_BANDING then
		-- decide which cysteine bands to enable and which to disable
		EnableAndDisableCysteineBands( )
	end
   
    -- decide whether to perform expansion
    if EXPAND_BEFORE_ACTION and randomBool( PROB_USE_EXPANDER ) then 
        PerformEveryXToEveryYExpansion()
        didExpand = true
    end

    -- non-band actions (rebuild, idealize, mutate)
    if algName == "SpotMutator" then
        changeSucceeded = PerformMutatorAction( idx )
    elseif algName == "IdealizePoint" then
        changeSucceeded = PerformIdealizePointAction( idx )
    elseif algName == "IdealizeInRange" then
        changeSucceeded = PerformIdealizeRangeAction( idx )
        
    elseif algName == "RebuildSS" then
        changeSucceeded = PerformRebuildSSAction( idx )
        isRebuild = true
    elseif algName == "RebuildInRange" then
        changeSucceeded = PerformRebuildRangeAction( idx )
        isRebuild = true
    elseif algName == "RebuildSimplePoint" then
        changeSucceeded = PerformRebuildSimplePointAction( idx )
        isRebuild = true
    elseif algName == "RebuildSimpleRange" then
        changeSucceeded = PerformRebuildSimpleRangeAction( idx )
        isRebuild = true
    elseif algName == "RebuildSimple2Point" then
        changeSucceeded = PerformRebuildSimple2PointAction( idx )
        isRebuild = true
    elseif algName == "RebuildSimple2Range" then
        changeSucceeded = PerformRebuildSimple2RangeAction( idx )
        isRebuild = true
    elseif algName == "RebuildSimple3Range" then
        changeSucceeded = PerformRebuildSimple3RangeAction( idx ) 
        isRebuild = true
    elseif algName == "RebuildFromQuaker" then
        changeSucceeded = PerformRebuildFromQuakerAction( idx ) 
        isRebuild = true
    elseif algName == "RebuildBitspawn" then
        changeSucceeded = PerformRebuildBitspawnAction( idx ) 
        isRebuild = true
    elseif algName == "RebuildBitspawnLowCI" then
        changeSucceeded = PerformRebuildBitspawnLowCI( idx ) 
        isRebuild = true

    -- simple-task banders
    elseif algName == "BIS" then
        changeSucceeded = PerformSingleBandInSpaceAction( idx )
    elseif algName == "BandBetween" then
        changeSucceeded = PerformBandToRandomSegAction( idx )
    elseif algName == "RotatorBands" then
        changeSucceeded = PerformRotatorAction( idx )
    elseif algName == "PusherBands" then
        changeSucceeded = PerformPusherAction( idx )
    elseif algName == "MoverBands" then
        changeSucceeded = PerformMoverAction( idx )
    elseif algName == "FlexerBands" then
        changeSucceeded = PerformFlexerAction( idx )
    elseif algName == "SheetStraightener" then
        changeSucceeded = PerformSheetStraightenerAction( idx )
    elseif algName == "HelixUniform" then
        changeSucceeded = PerformHelixUniformAction( idx )
    elseif algName == "HelixCompressor" then
        changeSucceeded = PerformHelixCompressorAction( idx )
    elseif algName == "HelixFixer" then
        changeSucceeded = PerformHelixFixerAction( idx )
    elseif algName == "SheetSewer" then
        changeSucceeded = PerformSheetSewerAction( idx )
    elseif algName == "GeneralHinge" then
        changeSucceeded = PerformGeneralHingeAction( idx )
    elseif algName == "GlycineHinge" then
        changeSucceeded = PerformGlycineHingeAction( idx )
    elseif algName == "TailWhip" then
        changeSucceeded = PerformTailWhipAction( idx )
    elseif algName == "PushPull" then
        changeSucceeded = PerformPushPullAction( idx )
    elseif algName == "VoidCrusher" then
        changeSucceeded = PerformVoidCrusherAction( idx )
        
    -- large-scale banders
    elseif algName == "MultiBander" then
        changeSucceeded = PerformMultiBanderAction( ) 
    elseif algName == "CentralBander" then
        changeSucceeded = PerformCentralBanderAction( ) 
    elseif algName =="EveryXToY" then
        changeSucceeded = PerformEveryXToEveryYBandAction( )
    elseif algName == "GentleComp" then
        changeSucceeded = PerformGentleCompAction( )
    elseif algName == "GentleCompLowCI" then
        changeSucceeded = PerformGentleCompLowCIAction( )
    elseif algName == "CoreStirrer" then
        changeSucceeded = PerformCoreStirrerAction( )
    elseif algName == "LSQuaker" then
        changeSucceeded = PerformLSQuakeAction( ) 
    elseif algName == "LocalQuake" then
        changeSucceeded = PerformLocalQuakeAction( idx )

    -- special contact map actions
    elseif algName == "ContactMapper" then
        changeSucceeded = PerformContactMapperAction( )
    elseif algName == "ContactMapTwoRegions" then
        changeSucceeded = PerformCMTwoRegionBandAction( idx )
    elseif algName == "ContactTwoSegs" then
        changeSucceeded = PerformSingleContactMapAction( idx )
    elseif algName == "SomeContactMapBands" then
        changeSucceeded = PerformSomeContactMapBandsAction( )
    elseif algName == "ContactMapForSeg" then
        changeSucceeded = PerformContactMapBandsForSegAction( idx )
    elseif algName == "PushPullContactMap" then
        changeSucceeded = PerformPushPullContactMapAction( idx )
    elseif algName == "VoidCrusherContactMap" then
        changeSucceeded = PerformVoidCrusherContactMapAction( idx )
    elseif algName == "ContactMapRebuildRange" then
        changeSucceeded = PerformContactMapRebuildRangeAction( idx )
        isRebuild = true
    elseif algName == "ContactMapRebuildSS" then
        changeSucceeded = PerformContactMapRebuildSSAction( idx )
        isRebuild = true
    elseif algName == "ContactMapRebuildAll" then
        changeSucceeded = PerformContactMapRebuildAllAction( )
        isRebuild = true
        
    else
        print( "Error - invalid action: ".. algName )
        return false
    end

    -- cleanup on aisle 6
    ResetFrozenness( )
    DelBands( )
    ResetCuts( )
    selection.DeselectAll( )
    SetCI( 1.0 )

    -- A bit of record-keeping for failures
    if not changeSucceeded then 
        -- heavily immovable proteins can fail rebuilds forever. 
        -- We need to be able to bail, so increment global failure variable
        if isRebuild then 
            rebuild_fails = rebuild_fails + 1
            print( "CHANGE FAILED - RB(" .. rebuild_fails .. ")" ) 
        else
            print( "CHANGE FAILED - " .. algName )
        end
        
        return false
    end

    -- closeout wiggles and such.
    if didExpand then
        algName = "Exp" .. algName            -- TEMPORARY
        
        -- undo our extra action at the start of it all
        recentbest.Save( )
        ShakeOrMutate( SHAKESIDECHAIN_ITERS )
        recentbest.Restore( )
        LongWiggle( )
        recentbest.Restore( )
    end

    -- Always do a closeout Shake/Mutate + LongWiggle	
    recentbest.Save( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    recentbest.Restore( )
    LongWiggle( )
    recentbest.Restore( )    

    -- if a fuse is going to happen, now is the time
    if USE_FUSEWIGGLE then 
        DoFuseWiggle( )
    end

    return changeSucceeded
end


----------------------------------------------------------------
-- MULTITRY FUNCTIONS
----------------------------------------------------------------
function SaveTopBests( TopBestSlots )
    -- todo: consider demanding a pose be significantly different 
    -- before adding (if not, but better score, then just replace)
    save.Quicksave( QS_Swapslot )
    currScore = getScore()
    for i = 1, #TopBestSlots do
        save.Quickload( TopBestSlots[ i ] )
        savScore = getScore()
        if currScore > savScore then
            for j=#TopBestSlots-1, i, -1 do
                save.Quickload( TopBestSlots[j] )
                save.Quicksave( TopBestSlots[j+1] )
            end
            save.Quickload( QS_Swapslot )
            save.Quicksave( TopBestSlots[ i ] ) -- save new pose
            break
        end
    end
    save.Quickload( QS_Swapslot )
end
function GetBestFromList( TopBests, wantReport )
    -- it's already sorted, so we'll load the first one
    if wantReport then
        for i=1, #TopBests do
            save.Quickload( TopBests[ i ] )
            print( "Score #" .. i .. " = " .. getScore() )
        end
    end
    save.Quickload( TopBests[1] )
end
function SetupTopBests( TopBests, NumberOfSlotsWanted )
    if NumberOfSlotsWanted > QS_MultiStartSaveLast - QS_MultiStartSave1 + 1 then
        print( "!!! Asked for more slots than allowed max" )
        return false
    end
    for i=1, NumberOfSlotsWanted do
        TopBests[#TopBests+1] = QS_MultiStartSave1 + i - 1
        save.Quicksave( QS_MultiStartSave1 + i - 1)
    end
    return true
end

function PerformMultiTryOperation( 
            NumberOfStarts, NumberOfActionsPerTry,
            startIdx, endIdx, ActionFunction, abbrBase )
    if abbrBase == nil then abbrBase = "" end
    
    print( "PerformMultiTryOperation for "..NumberOfStarts.." starts with "..
           NumberOfActionsPerTry.." actions per try" )
    BestScore=getScore( )
    LoadBest( )
    bests = {}
    TopBests = {}
    
    if not SetupTopBests( TopBests, NumberOfStarts ) then return false end

    save.Quicksave( QS_TopLevelTempMT )
    save.Quicksave( QS_TopLevelTempM2 )

    for k=1, NumberOfStarts do
        CurrentStart = k
        save.Quickload( QS_TopLevelTempMT )   -- always start over from initial state
        SaveBest( true, QS_TopLevelTempM2 )   -- best so far
        --print( "Starting try " .. k.. " of " .. NumberOfStarts ..  ": Score="..TrimNum(getScore( )))

        local initscore = getScore( )
        for j=1, NumberOfActionsPerTry do
            if NumberOfActionsPerTry > 1 then
                abbr = abbrBase .. "-MT" .. NumberOfStarts .. "." .. NumberOfActionsPerTry.. ":" .. k .. "." .. j
            else
                abbr = abbrBase .. "-MT" .. NumberOfStarts .. ":" .. k
            end
            save.Quickload( QS_TopLevelTempM2 )
            local currScore = getScore( )

            ActionFunction( startIdx, endIdx, QS_TopLevelTempM2, abbr )
            -- I expect ActionFunction to store its best result in QS..M2,
            -- but I do not expect it to load that result for me.
            save.Quickload( QS_TopLevelTempM2 ) 
            local endScore = getScore( )
            if endScore > currScore then
                save.Quicksave( QS_TopLevelTempM2 )
                currScore = endScore
            end
            if endScore > BestScore then
                BestScore = endScore
            end
            SaveBest( false, QS_TopLevelTempM2 )
            SaveBest( )
        end    -- FOR j

        save.Quickload( QS_TopLevelTempM2 )
        print( abbr .. ": Finishing - Score: "..getScore( )  )
        SaveTopBests( TopBests )
    end       -- FOR k

    GetBestFromList( TopBests, true )
    print( "Taking result with score " .. TrimNum( getScore() ))
    return true
end



----------------------------------------------------------------------
-- SCRIPT-SPECIFIC FUNCTIONS
----------------------------------------------------------------------

----------------------------------------------------------------------
-- ROUTINES FOR CHOOSING ACTION COLLECTIONS
----------------------------------------------------------------------
function AddRegionBanderActionsToList(allowedActions, idx )
    allowedActions[ #allowedActions + 1 ] = "PusherBands"
    allowedActions[ #allowedActions + 1 ] = "MoverBands"
    
    local startSeg, endSeg = GetRegionForOperation( idx )
    local regionLength = endSeg - startSeg + 1
    local isSheetRegion = structure.GetSecondaryStructure( idx ) == 'E'
    local isHelixRegion = structure.GetSecondaryStructure( idx ) == 'H' or isPossibleHelix( idx )
    
    -- Decide what other actions make sense to attempt 
    if (not isHelixRegion) and idx >= 5 and idx <= segCt - 5 then
        -- GeneralHinge tends to be useless
        allowedActions[ #allowedActions + 1 ] = "GeneralHinge"
    end
    if (idx > 3) and (idx < segCt-3) and (structure.GetAminoAcid( idx ) == 'g') then
         allowedActions[ #allowedActions + 1 ] = "GlycineHinge"
    end
    if isHelixRegion then
        -- rotator is just another "mover" if aas are not coiled as in a helix
        allowedActions[ #allowedActions + 1 ] = "RotatorBands"
        if regionLength > 6 then 
            -- actions that only make sense for helices 
            --    (not just "curls" of helix-ness)
            allowedActions[ #allowedActions + 1] = "HelixUniform"
            allowedActions[ #allowedActions + 1] = "HelixCompressor" 
            allowedActions[ #allowedActions + 1] = "HelixFixer"
        end
    end
    if isSheetRegion and #allSheets > 1 then
        -- SheetStraightener never turns out well...
        -- allowedActions[ #allowedActions + 1 ] = "SheetStraightener"
        allowedActions[ #allowedActions + 1] = "SheetSewer"
    end
    if regionLength >= 5 and idx >= 4 and idx <= segCt - 4 then
        allowedActions[ #allowedActions + 1] = "FlexerBands"
    end

    return
end


function PickMidgameAction( idx )
    if randomBool( PROB_TARGETED ) then
        idx = randomLowScoringSeg( random( segCt/5 ))
    else
        if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    end
    
    if idx == nil then idx = randomMovableSeg( ) end
    local allowedActions = { 
        "PushPull", 
        "VoidCrusher", 
        "GentleComp"
    }
    if (idx > 3) and (idx < segCt-3) and (structure.GetAminoAcid( idx ) == 'g') then
         allowedActions[ #allowedActions + 1 ] = "GlycineHinge"
    end
    if #NonLockedSegList > 0 then
        allowedActions[ #allowedActions + 1] = "RebuildInRange"
        allowedActions[ #allowedActions + 1] = "RebuildFromQuaker"
        allowedActions[ #allowedActions + 1] = "IdealizeInRange"
    end 
    AddRegionBanderActionsToList( allowedActions, idx )

    if PuzzleHasContactMap then 
        allowedActions[ #allowedActions + 1] = "ContactMapper" 
        allowedActions[ #allowedActions + 1] = "SomeContactMapBands" 
        allowedActions[ #allowedActions + 1] = "ContactMapTwoRegions"
        allowedActions[ #allowedActions + 1] = "ContactMapRebuildRange"
        allowedActions[ #allowedActions + 1] = "ContactMapRebuildSS"
        allowedActions[ #allowedActions + 1] = "ContactTwoSegs"
        if SegHasContactData( idx, CONTACTMAP_THRESHOLD ) then
            allowedActions[ #allowedActions + 1] = "ContactMapForSeg"
        end
    end

    return PickActionFromList( allowedActions, idx )
end

function PickEarlyAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    local allowedActions = { }
        
    allowedActions[ #allowedActions + 1] = "GentleComp"
    allowedActions[ #allowedActions + 1] = "PushPull"
    allowedActions[ #allowedActions + 1] = "GentleCompLowCI"
    allowedActions[ #allowedActions + 1] = "VoidCrusher"
    
    if #NonLockedSegList > 0 then
        allowedActions[ #allowedActions + 1] = "RebuildSimpleRange"
        allowedActions[ #allowedActions + 1] = "RebuildBitspawn"
        allowedActions[ #allowedActions + 1] = "RebuildBitspawnLowCI"
    end
    AddRegionBanderActionsToList( allowedActions, idx ) -- push, move, flex, helix, sheet
    
    if PuzzleHasContactMap then
        allowedActions[ #allowedActions + 1] = "ContactMapper"
        allowedActions[ #allowedActions + 1] = "SomeContactMapBands"
        allowedActions[ #allowedActions + 1] = "ContactMapTwoRegions"

        if #NonLockedSegList > 0 then
            allowedActions[ #allowedActions + 1] = "ContactMapRebuildRange"
        end

        if SegHasContactData( idx, CONTACTMAP_THRESHOLD ) then
            allowedActions[ #allowedActions + 1] = "PushPullContactMap"
            allowedActions[ #allowedActions + 1] = "VoidCrusherContactMap"
            allowedActions[ #allowedActions + 1] = "ContactMapForSeg"
        end
    end
    return PickActionFromList( allowedActions, idx )
end

function PickRebuildAction( idx )
    changeSucceeded = false
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    local allowedActions = { 
        "RebuildInRange", 
        "RebuildSimpleRange", 
        "RebuildSimple2Range", 
        "RebuildSimple3Range", 
        "RebuildFromQuaker", 
        "RebuildBitspawn", 
        "RebuildBitspawnLowCI", 
        "RebuildSS"
    }

    return PickActionFromList( allowedActions, idx )
end

function PickLowCIAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    local allowedActions = { }
        
    allowedActions[ #allowedActions + 1] = "GentleCompLowCI"
    allowedActions[ #allowedActions + 1] = "VoidCrusher"
    
    if #NonLockedSegList > 0 then
        allowedActions[ #allowedActions + 1] = "RebuildSimpleRange"
        allowedActions[ #allowedActions + 1] = "RebuildBitspawnLowCI"
    end
    AddRegionBanderActionsToList( allowedActions, idx ) -- push, move, flex, helix, sheet
    
    if PuzzleHasContactMap then
        allowedActions[ #allowedActions + 1] = "ContactMapper"
        allowedActions[ #allowedActions + 1] = "SomeContactMapBands"
        allowedActions[ #allowedActions + 1] = "ContactMapTwoRegions"

        if #NonLockedSegList > 0 then
            allowedActions[ #allowedActions + 1] = "ContactMapRebuildRange"
        end

        if SegHasContactData( idx, CONTACTMAP_THRESHOLD ) then
            allowedActions[ #allowedActions + 1] = "PushPullContactMap"
            allowedActions[ #allowedActions + 1] = "VoidCrusherContactMap"
            allowedActions[ #allowedActions + 1] = "ContactMapForSeg"
        end
    end
    return PickActionFromList( allowedActions, idx )
end

function PickLateAction( idx )
    if idx == nil or idx == 0 then idx = randomMovableSeg( ) end
    local allowedActions = { 
        "CentralBander", 
        "EveryXToY", 
        "PushPull", 
        "VoidCrusher", 
        "BIS", 
        "BandBetween",
        "IdealizePoint", 
    }
    if #NonLockedSegList > 0 then
        allowedActions[ #allowedActions + 1] = "RebuildInRange"
        allowedActions[ #allowedActions + 1] = "IdealizeInRange"
        allowedActions[ #allowedActions + 1] = "RebuildSimplePoint"
        allowedActions[ #allowedActions + 1] = "RebuildSimple2Point"
    end
    AddRegionBanderActionsToList( allowedActions, idx )
    
    if PuzzleHasContactMap and SegHasContactData( idx, CONTACTMAP_THRESHOLD ) then
        allowedActions[ #allowedActions + 1] = "ContactTwoSegs"
        allowedActions[ #allowedActions + 1] = "ContactMapTwoRegions"
        allowedActions[ #allowedActions + 1] = "ContactMapRebuildRange"
        allowedActions[ #allowedActions + 1] = "ContactMapRebuildSS"
        if SegHasContactData( idx, CONTACTMAP_THRESHOLD ) then
            allowedActions[ #allowedActions + 1] = "ContactMapForSeg"
        end
    end
    
    return PickActionFromList( allowedActions, idx )
end

function PickContactMapAction( idx )
    if idx == nil or idx == 0 then idx = randomContactableSeg( ) end
    local allowedActions = { 
        "ContactMapper", 
        "SomeContactMapBands", 
        "ContactMapForSeg", 
        "ContactTwoSegs", 
        "ContactMapTwoRegions",
        "ContactMapRebuildRange",
        "ContactMapRebuildSS",
        "ContactMapRebuildAll",
    }
    return PickActionFromList( allowedActions, idx )
end


----------------------------------------------------------------
-- SA's Neighbor generator
----------------------------------------------------------------
function CreateNeighborPose( threshold, NeighborGenerationFunction )
    local oldCreditBest = creditbest.GetScore( )  -- to find out if we've missed transients
    local oldScoreMat = {}
    local oldQuickScore = getQuickScore( )
    local oldScore = getScore( )
    ConstructScoreMatrix( oldScoreMat )

    SetCI( 1.0 )
    DelBands( )
    selection.DeselectAll( )
    ResetFrozenness( )

    recentbest.Save( )
    local rbscore = getRBScore( )

    save.Quicksave( QS_NeighborTempSA )
    local goodNeighbor = false
    while not goodNeighbor do
        local oldTime = os.time( )
        local changeSucceeded = false
        goodNeighbor = false

        save.Quickload( QS_NeighborTempSA )
        changeSucceeded = NeighborGenerationFunction( 0 )

        if changeSucceeded then
            local newScore = getScore( )
            goodNeighbor = IsProteinChanged( oldQuickScore, oldScoreMat, threshold )
            if not goodNeighbor then
                AddFailToStats( )
                print( "NEIGHBOR REJECTED (score:"..newScore..")")
            else
                local deltaScore = newScore - oldScore
                local deltaTime = os.time() - oldTime
                AddActionToStats( deltaScore, deltaTime )
            end
        else
            AddFailToStats( )
        end
        if rebuild_fails > MAX_REBUILD_FAILS then
            -- sometimes there are no successful ways to do a rebuild, we must bail
            rebuild_fails = 0
            return false
        end
    end   -- END while not goodNeighbor

    -- Look over the result
    local newScore = getScore( )
    local newCreditBest = creditbest.GetScore( )
    if newCreditBest ~= oldCreditBest and newCreditBest ~= newScore then
        -- we missed a transient
        --    trying to regain it can be problematic for some wigglers 
        --    (eg WiggleZlbNicely) because the creditbest pose may have bands.
        print( "!!! CreditBest changed: old="..oldCreditBest.." new="..newCreditBest )
    end

    DebugPrint( "Neighbor Score: " .. newScore)
    return true
end

function DoFuseWiggle( )
    oldScore = getQuickScore( )
    oldTime = os.time( )
    local r = random( ) 
    if r < PROB_FUSEWIGGLE_LOCAL then
        print( "IterativeWiggleLocalByChunk" )
        algName = "IterativeWiggleLocalByChunk"
        IterativeWiggleLocalByChunk( 1, segCt, FUSEWIGGLE_GAIN, random( 3, 6 ), false )
    elseif r < PROB_FUSEWIGGLE_LOCAL + PROB_FUSEWIGGLE_BLUEFUSE then
        print( "BlueFuse" )
        algName = "BlueFuse"
        BlueFuse( )
    else   -- PROB_FUSEWIGGLE_CUTFUSE
        print( "IterativeCutAndWiggle")
        algName = "IterativeCutAndWiggle"
        IterativeCutAndWiggle( FUSEWIGGLE_ITERS, FUSEWIGGLE_GAIN, random( 3,11 ) )
        for i=1, segCt do   -- deliberate non-use of ResetCuts
            DeleteCutAtIndex( i )
        end
    end
    local deltaScore = getQuickScore( ) - oldScore
    local deltaTime = os.time( ) - oldTime
    AddActionToStats( deltaScore, deltaTime )
end


----------------------------------------------------------------------
--
---- SIMULATED ANNEALING ALGORITHM FUNCTIONS
--
----------------------------------------------------------------------
-- Real research could be done (rosetta folks have done) on what a "good"
-- Temperature function is for these purposes. 
-- 1/x does well enough for early functions. 
-- Later I prefer Middling or Zero.
-- Note that "num" goes from 0 to 1 over the life of the SA function.
function T_OneOverX( num )
    return ( 1.0 / (num + 0.01) ) - ( 1.0 / 1.01 )
end
function T_OneOverX2( num ) 
    x = ( 1.0 / (num + 0.001) ) - 1.0
    return x * x
end
function T_MiddlingX( num ) 
    return 10.0 * (1.0 - num)
end
function T_LowX( num ) 
    return 1.0 - num
end
function T_ColdX( num ) 
    return 0.1 * (1.0 - num)
end
function T_Zero( num ) 
    return 0.0  -- causes SA to drop to normal march-forward behavior.
end
-- this is the standard probability function in Simulated Annealing algorithms
function P( ecurr, enew, T )  
    if enew < ecurr then return 1.0 end    -- always accept increased scores
    if T == 0.0 then return 0.0 end        -- always reject decreased scores if temp = 0
    return math.exp( (ecurr - enew)/T )
end
-- must negate score to function as a Simulated Annealing energy function
function E( )  
    return 8000.0 - getScore( )    -- 8000 is irrelevant, just for scaling
end


function SA( kmax, emax, failNoGainPt, NeighborGenerationFunction, 
             TempFunc, qs_best, abbr )  
    -- emax is expected to "never happen"
    if TempFunc == nil then TempFunc = T_OneOverX end
    local estart = E( )
    local ecurr = estart
    save.Quicksave( QS_TopLevelTempSA )  -- a private slot
    local ebest = ecurr 
    local startTime = os.time( )

    local k = 0.0 
    while k < kmax and ecurr > emax do
        save.Quickload( QS_TopLevelTempSA )

        local T = TempFunc( k / kmax )
        print("#" ..abbr .. k+1 ..
              "  tp=" .. TrimNum( T )..
              " sc=" .. TrimNum( getScore( ) )..
              " time=" .. os.time( )-startTime )
        local succeeded = CreateNeighborPose( 0.1, NeighborGenerationFunction )
        if not succeeded then return k end
        local enew = E( )
        
        local p = P(ecurr, enew, T)
        local r = random( )
        if p > r or enew < ecurr then
            if enew < ecurr then
                print( "++++ accepting state: score " ..TrimNum( getScore()))
                if enew < ebest then    
                    ebest = enew
                    SaveBest( false, qs_best )
                end
            elseif enew > ecurr then
                print( "++++ overriding state: score " ..TrimNum( getScore()))
            else
                print( "++++ same-score state: score " ..TrimNum( getScore()))
            end
            save.Quicksave( QS_TopLevelTempSA )
            ecurr = enew 
        end
        if k + 1 >= failNoGainPt and k + 1 < kmax and ebest >= estart then
            print( "---- Reached step "..(k+1).." with no improvement - bailing" )
            LoadBest( qs_best )
            return k+1    -- give up: we've not achieved any improvement and we think we never will
        end
        k = k + 1.0 
    end
    LoadBest( qs_best )
    return k
end

function PerformSA( failNoGainPt, ctRestarts, 
                    stepsPerRun, minGainPerRun, 
                    NeighborGenerationFunction, 
                    TempFunc, quitIfTooManyFails, 
                    qs_best, abbr )
    if minGainPerRun == nil then minGainPerRun = 0 end
    local saInitialTime = os.time( )
    local saInitialScore = getScore( )
    local bestScore = saInitialScore
    local newscore = 0.0
    local stepCount = 0
    SaveBest( false, qs_best )
    if abbr == nil then 
        abbr = "SA" 
    else
        abbr = abbr .. "-SA"
    end
    if ctRestarts > 1 then
        abbr = abbr .. ctRestarts .. "."
    end
    abbr = abbr .. stepsPerRun .. ":"
    
    local ctFailedRuns = 0
    local failPt = stepsPerRun
    if quitIfTooManyFails then 
        failPt = failNoGainPt 
    else 
        failPt = stepsPerRun 
    end

    for i=1, ctRestarts do
        LoadBest( qs_best )
        local saRunStartScore = getScore( )
        local saRunStartTime = os.time( )
        print( "START: SA  score=" ..TrimNum( saRunStartScore )..
               "  time=" ..saRunStartTime - saInitialTime )
        
        abbrOut = abbr
        if ctRestarts > 1 then
            abbrOut = abbrOut .. i .. "."
        end

        stepCount = stepCount + SA( stepsPerRun, -999999.0, failPt, 
                                    NeighborGenerationFunction, TempFunc, 
                                    qs_best, abbrOut )

        LoadBest( qs_best )
        newscore = getScore( )
        gain = newscore - bestScore
        if newscore > bestScore then
            bestScore = newscore
        end
        print( "FINISH: SA score=" ..TrimNum( bestScore )..
               " gain=" ..TrimNum( gain )..
               " time=" ..os.time( ) - saRunStartTime )
        
        if quitIfTooManyFails and gain < minGainPerRun then
            -- this run wasn't profitable...
            ctFailedRuns = ctFailedRuns + 1
            if ctFailedRuns > 3 and i < ctFailedRuns + 2 then
                -- ok, it was fun while it lasted, but this is getting silly
                print( "Too many failed runs - bailing" )
                break
            end
        end
    end
    LoadBest( qs_best )
    if ctRestarts > 1 then
        print( "FINAL SA: steps=" ..stepCount..
               "  gain=" ..TrimNum( getScore( ) - saInitialScore)..
               "  time=" ..os.time( ) - saInitialTime )
    end
    return stepCount
end


----------------------------------------------------------------
-- SA OPERATION PERFORMING FUNCTIONS
----------------------------------------------------------------
function PerformEarlyOperation()
    -- no expander at first, some later. never fuses, always short runs
    --   no expander at first is to go quicker.
    EXPAND_BEFORE_ACTION = false 
    USE_FUSEWIGGLE = false  
    stepsEarly = 0

    -- EARLY  without expander
    EXPAND_BEFORE_ACTION = false
    tempFunc = T_OneOverX 
    runcount = SA_NUMSTARTS_EARLY   
    stepcount = SA_STEPSPERRUN_EARLY
    print( "SA-Early for " .. stepcount ..
           " steps and " .. runcount ..
           " runs; do not fail out" )
           
    stepsEarly = PerformSA( stepcount, 
               runcount, stepcount, SA_MINGAIN, 
               PickEarlyAction, tempFunc, false, 
               QS_Best, "EA1" )
   
    --  EARLY   with expander
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION 
    PROB_USE_EXPANDER = PROB_USE_EXPANDER_EARLY2 
    tempFunc = T_MiddlingX          -- T_OneOverX T_MiddlingX
    runcount = SA_NUMSTARTS_EARLY2   
    stepcount = SA_STEPSPERRUN_EARLY2
    print( "SA-Early for " .. stepcount ..
           " steps and " .. runcount ..
           " runs; do not fail out" )
           
    stepsEarly = stepsEarly + 
                PerformSA( stepcount, 
                           runcount, stepcount, SA_MINGAIN, 
                           PickEarlyAction, tempFunc, false, 
                           QS_Best, "EA2" )
    return stepsEarly
end

function PerformRebuildOperation( )
    -- no fuses, some expander, short runs, 2nd cool
    if #NonLockedSegList == 0 then
        print( "Cannot do rebuild with no movable segments" )
        return 0
    end

    EXPAND_BEFORE_ACTION = false
    USE_FUSEWIGGLE = false  
    stepsRebuild = 0

    -- REBUILD with expander
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION 
    PROB_USE_EXPANDER = PROB_USE_EXPANDER_REBUILD 
    tempFunc = T_OneOverX 
    runcount = SA_NUMSTARTS_REBUILD
    numSteps = SA_STEPSPERRUN_REBUILD 
    print( "SARebuild-1/x for " .. numSteps ..
           " steps and " .. runcount ..
           " runs; do not fail out" )
    stepsRebuild = PerformSA( numSteps, 
                runcount, numSteps, SA_MINGAIN, 
                PickRebuildAction, tempFunc, false, 
                QS_Best, "RB1" )
   

    -- REBUILD with expander and a cooler temperature
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION 
    PROB_USE_EXPANDER = PROB_USE_EXPANDER_REBUILD 
    tempFunc = T_MiddlingX 
    runcount = SA_NUMSTARTS_REBUILD2
    numSteps = SA_STEPSPERRUN_REBUILD2
    print( "SARebuild-Middling for " .. numSteps ..
           " steps and " .. runcount ..
           " runs; do not fail out" )
    stepsRebuild = stepsRebuild + PerformSA( numSteps, 
                runcount, numSteps, SA_MINGAIN, 
                PickRebuildAction, tempFunc, false, 
                QS_Best, "RB2" )
    return stepsRebuild
end

function PerformMidgameOperation()
    -- no fuses, yes expander, be grateful for anything we get.
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION
    
    PROB_USE_EXPANDER = PROB_USE_EXPANDER_MID1 
    print( "SA-Midgame1 for " .. SA_STEPSPERRUN_MID ..
           " steps and " .. SA_NUMSTARTS_MID ..
           " runs; do not fail out" )
           
    stepCt = PerformSA( SA_STEPSPERRUN_MID, 
               SA_NUMSTARTS_MID, SA_STEPSPERRUN_MID, SA_MINGAIN, 
               PickMidgameAction, T_Zero, false, 
               QS_Best, "MID1" )


    PROB_USE_EXPANDER = PROB_USE_EXPANDER_MID2   
    print( "SA-Midgame2 for " .. SA_STEPSPERRUN_MID2 ..
           " steps and " .. SA_NUMSTARTS_MID2 ..
           " runs; do not fail out" )
           
    stepCt = stepCt + 
             PerformSA( SA_STEPSPERRUN_MID2, 
               SA_NUMSTARTS_MID2, SA_STEPSPERRUN_MID2, SA_MINGAIN, 
               PickMidgameAction, T_Zero, false, 
               QS_Best, "MID2" )
    return stepCt
end

function PerformLateOperation( )
    -- yes fuses, yes expander, be grateful for anything we get.
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION
    PROB_USE_EXPANDER = PROB_USE_EXPANDER_LATE 
    print( "SA-Late for " .. SA_STEPSPERRUN_LATE ..
           " steps and " .. SA_NUMSTARTS_LATE ..
           " runs; do not fail out" )
    local qwg = QUICKWIGGLE_GAIN
    local swg = FUSEWIGGLE_GAIN
    QUICKWIGGLE_GAIN = 0.1
    FUSEWIGGLE_GAIN = 0.01
    stepCt = PerformSA( SA_STEPSPERRUN_LATE, 
                   SA_NUMSTARTS_LATE, SA_STEPSPERRUN_LATE, SA_MINGAIN, 
                   PickLateAction, T_Zero, false, 
                   QS_Best, "LATE"..SA_STEPSPERRUN_LATE )
    QUICKWIGGLE_GAIN = qwg 
    FUSEWIGGLE_GAIN = swg
    
    return stepCt
end

function PerformLowCIOperation( )
    -- Does not seem worth the time it takes for "stable" poses
    -- pre-expander, no fuses, CI<=0.5
    savCI = GetCI( )
    EXPAND_BEFORE_ACTION = ALLOW_EXPAND_FOR_ACTION
    USE_LOWCI_WIGGLES = true
    ExactSetCI( 0.20 )
    SCALE_MAXCI = 0.20
    print( "SA-LowCI20 for " .. SA_STEPSPERRUN_LOW20 ..
           " steps and " .. SA_NUMSTARTS_LOW20 ..
           " runs; fail if no gain by " .. SA_FAILIFNOGAIN )
    stepCt = PerformSA( SA_FAILIFNOGAIN, 
               SA_NUMSTARTS_LOW20, SA_STEPSPERRUN_LOW20, SA_MINGAIN, 
               PickLowCIAction, T_OneOverX, true, 
               QS_Best, "20CI"..SA_STEPSPERRUN_LOW20 )

    ExactSetCI( 0.50 )
    SCALE_MAXCI = 0.50
    print( "SA-LowCI50 for " .. SA_STEPSPERRUN_LOW50 ..
           " steps and " .. SA_NUMSTARTS_LOW50 ..
           " runs; fail if no gain by " .. SA_FAILIFNOGAIN )
    stepCt = stepCt + PerformSA( SA_FAILIFNOGAIN, 
               SA_NUMSTARTS_LOW50, SA_STEPSPERRUN_LOW50, SA_MINGAIN, 
               PickLowCIAction, T_OneOverX, true, 
               QS_Best, "50CI"..SA_STEPSPERRUN_LOW50 )
    USE_LOWCI_WIGGLES = false
    ExactSetCI( savCI )
    SCALE_MAXCI = savCI
    
    return stepCt
end

function PerformContactMapOperation( )
    -- not in basic set, must have been specifically asked for
    print( "SA-Contact for " .. SA_STEPSPERRUN_CMAP ..
           " steps and " .. SA_NUMSTARTS_CMAP ..
           " runs; fail if no gain by " .. SA_FAILIFNOGAIN )
    stepCt = PerformSA( SA_FAILIFNOGAIN, 
               SA_NUMSTARTS_CMAP, SA_STEPSPERRUN_CMAP, SA_MINGAIN, 
               PickContactMapAction, T_OneOverX, true, 
               QS_Best, "CT"..SA_STEPSPERRUN_CMAP )
    
    return stepCt
end

function PerformOperation( OperationChosen )
    local moInitialScore = getScore( )
    local moInitialTime = os.time( )
    local numSteps = 0
    local failNoGain = 0
    local stepCt = 0
    rebuild_fails = 0 -- global that gets incremented when a rebuild fails 

    SaveBest( )
    if rebuild_fails >= MAX_REBUILD_FAILS then
        print( "rebuild fail count exceeded in PerformOperation")
        rebuild_fails = 0
        return
    end
    
    print( "Performing ".. OperationChosen )
    
    if OperationChosen == "SASequence" then
        PerformSASequence( )
    elseif OperationChosen == "SAEarly" then
        stepCt = PerformEarlyOperation( )
        print( "SA-Early finished" )
    elseif OperationChosen == "SARebuild" then
        stepCt = PerformRebuildOperation( )
        print( "SARebuild finished" )
    elseif OperationChosen == "SAMidgame" then
        stepCt = PerformMidgameOperation()
        print( "SA-Midgame finished" )

    elseif OperationChosen == "SALowCI" then
        stepCt = PerformLowCIOperation( )
        print( "SA-LowCI finished" )
    elseif OperationChosen == "SALate" then
        stepCt = PerformLateOperation( )
        print( "SA-Late finished" )
    elseif OperationChosen == "SAContact" then
        stepCt = PerformContactMapOperation( )
        print( "SA-Contact finished" )
        
    else
        print( "!! Unknown operation" )
    end

    LoadBest( )
    local endScore = getScore( )
    return stepCt
end

function PerformSASequence( )    
    SaveBest( )
    print( "Begin SA Sequence" )
    scoreInit = getScore( )
    timeInit = os.time( )

    -- INIT WSWSW
    timeStart = os.time()
	EnableAndDisableCysteineBands( )
    LongWiggle( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    LongWiggle( )
    ShakeOrMutate( SHAKESIDECHAIN_ITERS )
    LongWiggle( )
    scoreWswsw = getScore( )
    timeWswsw = os.time( ) - timeStart
	DelBands()
    SaveBest()

    -- EARLY
    timeStart = os.time()
    stepsEarly = PerformEarlyOperation( )
    timeEarly = os.time() - timeStart
    scoreEarly = getScore( )
    save.Quicksave( 4 )
    print( "SASEQ-Early finished - stored in qs 4" )


    -- REBUILDS
    timeStart = os.time()
    stepsRebuild = PerformRebuildOperation( )
    timeRebuild = os.time() - timeStart
    scoreRebuild = getScore( )
    save.Quicksave( 5 )
    print( "SARebuild finished - stored in qs 5" )    
    
    -- MIDGAME
    timeStart = os.time()
    stepsMid = PerformMidgameOperation( )
    timeMid = os.time( ) - timeStart
    scoreMid = getScore( )
    print( "SAMidgame finished" )

    
    -- FINAL   -- I question caller's wisdom, but do as they asked...
    if ALLOW_FUSE then
        save.Quicksave( 6 )
        print( "SAMidgame result stored in qs 6" )
        
        USE_FUSEWIGGLE = true
        IterativeCutAndWiggle( FUSEWIGGLE_ITERS, FUSEWIGGLE_GAIN, random( 3,11 ))
        timeStart = os.time()
        stepsLate = PerformLateOperation( )
        timeLate = os.time( ) - timeStart
        scoreLate = getScore( )
    end

    print( "End SA Sequence: gain="..getScore( ) - scoreInit..
           " time="..os.time( ) - timeInit )
    print( "  Stage    Time    Steps   Score" )
    print( "  Init:    " ..PadNumString(0,6)..          "  " ..PadNumString(0,5).. "  " .. TrimNum(scoreInit))
    print( "  WSWSW:   " ..PadNumString(timeWswsw,6)..  "  " ..PadNumString(0,5).. "  " .. TrimNum(scoreWswsw))
    print( "  Early:   " ..PadNumString(timeEarly,6)..  "  " ..PadNumString(stepsEarly,5).. "  " .. TrimNum(scoreEarly))
    print( "  Rebuild: " ..PadNumString(timeRebuild,6).."  " ..PadNumString(stepsRebuild,5).. "  " .. TrimNum(scoreRebuild))
    print( "  Midgame: " ..PadNumString(timeMid,6)..    "  " ..PadNumString(stepsMid,5).. "  " .. TrimNum(scoreMid))
    if ALLOW_FUSE then
        print( "  Late:    " ..PadNumString(timeLate,6).. "  " ..PadNumString(stepsLate,5).. "  " .. TrimNum(scoreLate))
    end
end



-----------------------------------------------------------------------------
-- DIALOG BOX SEQUENCE FROM THE NETHERWORLD
-----------------------------------------------------------------------------
function PutSlowWiggleOptionsInBox( dlg )
    dlg.sepSWSA = dialog.AddLabel( "----------------Fuse Wiggles------------------------" )
    dlg.swgain = dialog.AddSlider( "Fuse-WiggleGain", FUSEWIGGLE_GAIN, 0, 0.50, 2 )
    dlg.comment = dialog.AddLabel( "Fuse-wiggle values are 'probabilities' - scaled to sum to 1" )
    dlg.wloc = dialog.AddSlider( "LocalWiggler", PROB_FUSEWIGGLE_LOCAL, 0, 1, 2 ) -- slow, poor
    dlg.wbfuse= dialog.AddSlider( "BlueFuse", PROB_FUSEWIGGLE_BLUEFUSE, 0, 1, 2 )  -- less slow, moderately good
    dlg.wcwfuse= dialog.AddSlider( "C&WFuse", PROB_FUSEWIGGLE_CUTFUSE, 0, 1, 2 ) -- slow, moderately good
end
function SetSlowWiggleOptions( dlg )
    FUSEWIGGLE_GAIN = dlg.swgain.value
    local t1 = dlg.wloc.value
    local t2 = dlg.wbfuse.value
    local t3 = dlg.wcwfuse.value
    if t1 + t2 + t3 == 0 then
        t1 = .25  -- not a problem that these don't sum to 1: they'll be scaled
        t2 = .25
        t3 = .25
    end
    PROB_FUSEWIGGLE_LOCAL =    t1 / ( t1 + t2 + t3 )
    PROB_FUSEWIGGLE_BLUEFUSE = t2 / ( t1 + t2 + t3 )
    PROB_FUSEWIGGLE_CUTFUSE =  t3 / ( t1 + t2 + t3 )

    DebugPrint( "SlowGain="..FUSEWIGGLE_GAIN )
    DebugPrint( "Prob_slowwiggle-local="..PROB_FUSEWIGGLE_LOCAL )
    DebugPrint( "Prob_slowwiggle-bluefuse="..PROB_FUSEWIGGLE_BLUEFUSE )
    DebugPrint( "Prob_slowwiggle-cutfuse="..PROB_FUSEWIGGLE_CUTFUSE )
end

function PutBanderOptionsInBox( dlg )
    dlg.sepBander = dialog.AddLabel( "----------------Bander Options---------------------" )
    dlg.bdcpci = dialog.AddSlider( "CorePulling CI", COREPULLING_CI, 0.01, 1.0, 2 )
    dlg.bdregstren = dialog.AddSlider( "RegionBandStrength", BAND_STRENGTH_REGION, 0.50, 2.0, 1 )
    dlg.bdwigtargdelta = dialog.AddSlider( "BWig ScoreDelta", BANDEDWIGGLE_TARGET, 50, 1000, 0 )
    dlg.bdwigminexpand = dialog.AddSlider( "MinExpand", BAND_MINEXPAND, 1.00, 1.10, 2 )
    dlg.bdwigmaxexpand = dialog.AddSlider( "MaxExpand", BAND_MAXEXPAND, 1.10, 3.00, 2 )
    dlg.bdwigmincompr = dialog.AddSlider( "MinCompress", BAND_MINCOMPRESS, 0.85, 1.00, 2 )
    dlg.bdwigmaxcompr = dialog.AddSlider( "MaxCompress", BAND_MAXCOMPRESS, 0.25, 0.85, 2 )
end
function SetBanderOptions( dlg )
    COREPULLING_CI = dlg.bdcpci.value
    BAND_STRENGTH_REGION = dlg.bdregstren.value
    BANDEDWIGGLE_TARGET = dlg.bdwigtargdelta.value
    BAND_MINCOMPRESS = dlg.bdwigmincompr.value
    BAND_MAXCOMPRESS = dlg.bdwigmaxcompr.value
    BAND_MINEXPAND = dlg.bdwigminexpand.value
    BAND_MAXEXPAND = dlg.bdwigmaxexpand.value

    DebugPrint("Core-pulling CI="..COREPULLING_CI )
    DebugPrint("Region-bander strength="..BAND_STRENGTH_REGION )
    DebugPrint("Banded Wiggler Target Delta="..BANDEDWIGGLE_TARGET )
    DebugPrint("Max Compress="..BAND_MINCOMPRESS)
    DebugPrint("Min Compress="..BAND_MAXCOMPRESS)
    DebugPrint("Min Expand="..BAND_MINEXPAND )
    DebugPrint("Max Expand="..BAND_MAXEXPAND )
end

function PutIdealizeOptionsInBox( dlg )
    dlg.sepIdealizer = dialog.AddLabel( "--------------Idealizer Options--------------------" )
    dlg.idMaxRange = dialog.AddSlider( "MaxIdealizeRange", IDEALIZE_MAX_SEGLEN, 1, 20, 0 ) 
    dlg.idusecuts = dialog.AddCheckbox( "Allow cuts during idealize", IDEALIZE_USE_CUTS )
    dlg.idmusthealcuts = dialog.AddCheckbox( "Forbid unhealed cuts after idealize", IDEALIZE_FORBID_UNHEALED_CUTS )
end
function SetIdealizeOptions( dlg )
    IDEALIZE_MAX_SEGLEN = dlg.idMaxRange.value
    IDEALIZE_USE_CUTS = dlg.idusecuts.value
    IDEALIZE_FORBID_UNHEALED_CUTS = dlg.idmusthealcuts.value

    DebugPrint("Idealize max range="..IDEALIZE_MAX_SEGLEN )
    if IDEALIZE_USE_CUTS then DebugPrint("Cuts will be made") else DebugPrint("Cuts are forbidden") end
    if IDEALIZE_FORBID_UNHEALED_CUTS then DebugPrint("Cuts must heal to be accepted") else DebugPrint("Cuts permitted to survive") end
end

function PutGeneralRebuildOptionsInBox( dlg )
    dlg.sepRb = dialog.AddLabel( "---------------Rebuilder Options------------------" )
    dlg.rebMaxTries = dialog.AddSlider( "TriesPerStep", REBUILD_TRYCOUNT, 1, 20, 0 )
    dlg.rebFailsbeforequit = dialog.AddSlider( "FailsAllowed", REBUILD_FAILSBEFOREQUIT, 10, 100, 0 )
    dlg.rebuildMaxRange = dialog.AddSlider( "MaxRebuildRange", REBUILD_MAX_SEGLEN, 3, 20, 0 ) 
    dlg.rebWantIdealize = dialog.AddCheckbox( "Idealize", REBUILD_ADDIDEALIZE )
    dlg.rebIdealizeShouldUseCuts = dialog.AddCheckbox( "Idealize uses cuts", REBUILD_USECUTS ) 
end
function SetGeneralRebuildOptions( dlg )
    REBUILD_TRYCOUNT = dlg.rebMaxTries.value
    REBUILD_FAILSBEFOREQUIT = dlg.rebFailsbeforequit.value
    REBUILD_MAX_SEGLEN = dlg.rebuildMaxRange.value
    REBUILD_ADDIDEALIZE = dlg.rebWantIdealize.value
    REBUILD_USECUTS = dlg.rebIdealizeShouldUseCuts.value

    DebugPrint( "TriesPerStep="..REBUILD_TRYCOUNT )
    DebugPrint( "FailsBeforeQuit="..REBUILD_FAILSBEFOREQUIT )
    DebugPrint( "MaxRebuildSegmentLength="..REBUILD_MAX_SEGLEN )
    if REBUILD_ADDIDEALIZE then DebugPrint( "Rebuilds followed by idealize" ) else DebugPrint( "Rebuild do not idealize" ) end
    if REBUILD_USECUTS then DebugPrint( "Rebuild+idealize can use cuts" ) else DebugPrint( "Rebuild+idealize cannot use cuts" ) end
end
function PutLowCI20OptionsInBox( dlg )
    dlg.stepsperLCI20 = dialog.AddSlider( "Steps20CI", SA_STEPSPERRUN_LOW20, 10, 200, 0 )
    dlg.startsLCI20 = dialog.AddSlider(   "Runs20CI", SA_NUMSTARTS_LOW20, 1, 20, 0 ) 
end
function SetLowCI20Options( dlg )
    SA_STEPSPERRUN_LOW20 = dlg.stepsperLCI20.value
    SA_NUMSTARTS_LOW20 = dlg.startsLCI20.value
    DebugPrint("LowCI steps,starts="..SA_STEPSPERRUN_LOW20..","..SA_NUMSTARTS_LOW20 )
end
function PutLowCI50OptionsInBox( dlg )
    dlg.stepsperLCI50 = dialog.AddSlider( "Steps50CI", SA_STEPSPERRUN_LOW50, 10, 200, 0 )
    dlg.startsLCI50 = dialog.AddSlider(   "Runs50CI", SA_NUMSTARTS_LOW50, 1, 20, 0 ) 
end
function SetLowCI50Options( dlg )
    SA_STEPSPERRUN_LOW50 = dlg.stepsperLCI50.value
    SA_NUMSTARTS_LOW50 = dlg.startsLCI50.value
    DebugPrint("LowCI steps,starts="..SA_STEPSPERRUN_LOW50..","..SA_NUMSTARTS_LOW50 )
end
function PutEarlyOptionsInBox( dlg )
    dlg.stepsperEa = dialog.AddSlider( "StepsEarly1", SA_STEPSPERRUN_EARLY, 10, 50, 0 )
    dlg.startsEa = dialog.AddSlider(   "RunsEarly1", SA_NUMSTARTS_EARLY, 1, 300, 0 ) 
    dlg.stepsperEa2 = dialog.AddSlider( "StepsEarly2", SA_STEPSPERRUN_EARLY2, 10, 50, 0 )
    dlg.startsEa2 = dialog.AddSlider(   "RunsEarly2", SA_NUMSTARTS_EARLY2, 1, 300, 0 ) 
    dlg.probExpandEa2 = dialog.AddSlider("ProbExpanderE2", PROB_USE_EXPANDER_EARLY2, 0, 1.0, 2 ) 
end
function SetEarlyOptions( dlg )
    SA_STEPSPERRUN_EARLY = dlg.stepsperEa.value
    SA_NUMSTARTS_EARLY = dlg.startsEa.value
    SA_STEPSPERRUN_EARLY2 = dlg.stepsperEa2.value
    SA_NUMSTARTS_EARLY2 = dlg.startsEa2.value
    PROB_USE_EXPANDER_EARLY2 = dlg.probExpandEa2.value
    DebugPrint("Early steps,starts="..SA_STEPSPERRUN_EARLY..","..SA_NUMSTARTS_EARLY )
    DebugPrint("Early2 steps,starts="..SA_STEPSPERRUN_EARLY2..","..SA_NUMSTARTS_EARLY2 )
    DebugPrint("Prob Early Expander = " .. PROB_USE_EXPANDER_EARLY2)
end
function PutRebuildOptionsInBox( dlg )
    dlg.stepsperRb = dialog.AddSlider( "StepsRebuild1", SA_STEPSPERRUN_REBUILD, 10, 50, 0 )
    dlg.startsRb = dialog.AddSlider(   "RunsRebuild1", SA_NUMSTARTS_REBUILD, 1, 300, 0 ) 
    dlg.stepsperRb2 = dialog.AddSlider( "StepsRebuild2", SA_STEPSPERRUN_REBUILD2, 10, 50, 0 )
    dlg.startsRb2 = dialog.AddSlider(   "RunsRebuild2", SA_NUMSTARTS_REBUILD2, 1, 300, 0 ) 
    dlg.probExpandRb = dialog.AddSlider("ProbExpanderR", PROB_USE_EXPANDER_REBUILD, 0, 1.0, 2 ) 
end
function SetRebuildOptions( dlg )
    SA_STEPSPERRUN_REBUILD = dlg.stepsperRb.value
    SA_NUMSTARTS_REBUILD = dlg.startsRb.value
    SA_STEPSPERRUN_REBUILD2 = dlg.stepsperRb2.value
    SA_NUMSTARTS_REBUILD2 = dlg.startsRb2.value
    PROB_USE_EXPANDER_REBUILD = dlg.probExpandRb.value
    DebugPrint("Rebuild steps,starts="..SA_STEPSPERRUN_REBUILD..","..SA_NUMSTARTS_REBUILD )
    DebugPrint("Rebuild2 steps,starts="..SA_STEPSPERRUN_REBUILD2..","..SA_NUMSTARTS_REBUILD2 )
    DebugPrint("Prob Rebuild Expander = " .. PROB_USE_EXPANDER_REBUILD)
end
function PutMidOptionsInBox( dlg )
    dlg.stepsperMv = dialog.AddSlider( "StepsMid1", SA_STEPSPERRUN_MID, 100, 500, 0 )
    dlg.startsMv = dialog.AddSlider(   "RunsMid1", SA_NUMSTARTS_MID, 1, 5, 0 ) 
    dlg.probExpandMv = dialog.AddSlider("ProbExpanderM1", PROB_USE_EXPANDER_MID1, 0, 1.0, 2 ) 
    dlg.stepsperMid = dialog.AddSlider("StepsMid2", SA_STEPSPERRUN_MID2, 100, 500, 0 )
    dlg.startsMid = dialog.AddSlider(  "RunsMid2", SA_NUMSTARTS_MID2, 1, 5, 0 ) 
    dlg.probExpandMid = dialog.AddSlider("ProbExpanderM2", PROB_USE_EXPANDER_MID2, 0, 1.0, 2 ) 
    dlg.probTarg = dialog.AddSlider( "Prob Mid2 Targeted", PROB_TARGETED, 0, 1.0, 2 )
    end
function SetMidOptions( dlg )
    SA_STEPSPERRUN_MID2 = dlg.stepsperMid.value
    SA_NUMSTARTS_MID2 = dlg.startsMid.value
    PROB_USE_EXPANDER_MID2 = dlg.probExpandMid.value
    SA_STEPSPERRUN_MID = dlg.stepsperMv.value
    SA_NUMSTARTS_MID = dlg.startsMv.value
    PROB_USE_EXPANDER_MID1 = dlg.probExpandRb.value
    PROB_TARGETED = dlg.probTarg.value

    DebugPrint("Mid1 steps,starts="..SA_STEPSPERRUN_MID..","..SA_NUMSTARTS_MID )
    DebugPrint("Prob Mid1 Expander = " .. PROB_USE_EXPANDER_MID1)
    DebugPrint("Mid2 steps,starts="..SA_STEPSPERRUN_MID2..","..SA_NUMSTARTS_MID2 )
    DebugPrint("Prob Mid2 Expander = " .. PROB_USE_EXPANDER_MID2)
    DebugPrint("Prob Mid2 uses targeted operations = " .. PROB_TARGETED)
end
function PutLateOptionsInBox( dlg )
    dlg.stepsperLt = dialog.AddSlider( "StepsLate", SA_STEPSPERRUN_LATE, 10, 500, 0 )
    dlg.startsLt = dialog.AddSlider(   "RunsLate", SA_NUMSTARTS_LATE, 1, 10, 0 ) 
    dlg.probExpandLt = dialog.AddSlider("ProbExpander", PROB_USE_EXPANDER_LATE, 0, 1.0, 2 ) 
end
function SetLateOptions( dlg )
    if ALLOW_FUSE then
        SA_STEPSPERRUN_LATE = dlg.stepsperLt.value
        SA_NUMSTARTS_LATE = dlg.startsLt.value
        PROB_USE_EXPANDER_LATE = dlg.probExpandLt.value
        DebugPrint("Late steps,starts="..SA_STEPSPERRUN_LATE..","..SA_NUMSTARTS_LATE )
        DebugPrint("Prob Late Expander = " .. PROB_USE_EXPANDER_LATE)
    end
end
function PutContactMapOptionsInBox( dlg )
    dlg.stepsperCMap = dialog.AddSlider( "StepsCMap", SA_STEPSPERRUN_CMAP, 10, 200, 0 )
    dlg.startsCMap = dialog.AddSlider(   "RunsLowCMap", SA_NUMSTARTS_CMAP, 1, 20, 0 ) 
end
function SetContactMapOptions( dlg )
    SA_STEPSPERRUN_CMAP = dlg.stepsperCMap.value
    SA_NUMSTARTS_CMAP = dlg.startsCMap.value
    DebugPrint("CMap steps,starts="..SA_STEPSPERRUN_CMAP..","..SA_NUMSTARTS_CMAP )
end
function PutGainOptionsInBox( dlg )
    dlg.samingain = dialog.AddSlider(  "MinGainForRuns", SA_MINGAIN, 0, 5, 1 )
    dlg.quitnogain = dialog.AddSlider( "StopIfNoGain", SA_FAILIFNOGAIN, 1, 200, 0 )
end
function SetGainOptions( dlg )
    SA_MINGAIN = dlg.samingain.value
    SA_FAILIFNOGAIN = dlg.quitnogain.value
    DebugPrint("Step to stop at if no gain="..SA_FAILIFNOGAIN ) 
    DebugPrint("SA halt if gain < "..SA_MINGAIN)
end

function PutSASequenceOptionsInBox( dlg )
    PutEarlyOptionsInBox( dlg )
    PutRebuildOptionsInBox( dlg )
    PutMidOptionsInBox( dlg )
    dlg.wantFuses = dialog.AddCheckbox( "Add final fuse run", false )
    PutLateOptionsInBox( dlg )
    PutGainOptionsInBox( dlg )
end
function SetSASequenceOptions( dlg )
    SetEarlyOptions( dlg )
    SetRebuildOptions( dlg )
    SetMidOptions( dlg )
    ALLOW_FUSE = dlg.wantFuses.value
    SetLateOptions( dlg )
    SetGainOptions( dlg )
end

function AddBasicSpecialOptions( dlg, includeBander, allowFuseWiggle )
    dlg.cmtSO = dialog.AddLabel( "----- Advanced SA Parameters -----" )
    dlg.IdRAdvancedPars= dialog.AddCheckbox( "Idealize/Rebuild", false )
    if allowFuseWiggle then
        dlg.WigglerPars= dialog.AddCheckbox( "Fuser", false )
    end
    if includeBander then
        dlg.BanderPars= dialog.AddCheckbox( "Bander", false )
    end
    dlg.maxCI = dialog.AddSlider( "MaxCI", SCALE_MAXCI, 0, 1, 2 )
    dlg.ok = dialog.AddButton( "OK", 1 )
    dlg.cancel = dialog.AddButton( "Cancel", 0 )
end
function SetBasicSpecialOptions( dlg, includeBander, allowFuseWiggle )
    if  not dlg.IdRAdvancedPars.value and 
        not ((not forbidFuseWiggle) and dlg.WigglerPars.value) and
        not (includeBander and dlg.BanderPars.value)
    then
        return
    end
    
    SCALE_MAXCI = dlg.maxCI.value
    DebugPrint( "Using MaxCI="..SCALE_MAXCI )
    
    dlg2 = dialog.CreateDialog( "Special" )
    if dlg.IdRAdvancedPars.value then
        PutIdealizeOptionsInBox( dlg2 )
        PutGeneralRebuildOptionsInBox( dlg2 )
    end
    if allowFuseWiggle and dlg.WigglerPars.value then
        PutSlowWiggleOptionsInBox( dlg2 )
    end
    if includeBander and dlg.BanderPars.value then
        PutBanderOptionsInBox( dlg2 )
    end
    dlg2.ok = dialog.AddButton( "OK", 1 )
    dlg2.cancel = dialog.AddButton( "Cancel", 0 )
    if dialog.Show( dlg2 ) == 0 then return false end
    
    if dlg.IdRAdvancedPars.value then
        SetIdealizeOptions( dlg2 )
        SetGeneralRebuildOptions( dlg2 )
    end
    if allowFuseWiggle and dlg.WigglerPars.value then
        SetSlowWiggleOptions( dlg2 )
    end
    if includeBander and dlg.BanderPars.value then
        SetBanderOptions( dlg2 )
    end
end
function PutCysteineOptionsInBox( dlg )
    str = ""
    for i=1, #cysteineSegs do
        str = str..cysteineSegs[i]
        if i < #cysteineSegs then
            str = str.. ","
        end
    end
    dlg.cmtCy = dialog.AddLabel( "Cysteines: "..str )
    str = ""
    for i=1, #cysteinePairs do
        pair = cysteinePairs[i]
        str = str..pair[1]..","..pair[2]
        if i < #cysteinePairs then
            str = str .. "; "
        end
    end
    dlg.cmtCy2 = dialog.AddLabel( "Guessed pairs: "..str )
    dlg.cyPairsWanted = dialog.AddSlider( "Pair count", #cysteinePairs, 0, #cysteinePairs, 0 )
    dlg.cyPairings = dialog.AddTextbox( "Pairings", str)
	if InitialBandCount > 0 then
		dlg.cmtCy3 = dialog.AddLabel( "Cysteine-bander cannot honor user bands" )
	end
    dlg.ok = dialog.AddButton( "OK", 1 )
    dlg.cancel = dialog.AddButton( "Cancel", 0 )
end

function parseCysteineStr( ctPairsWanted, str )
    local cysList = {}
    -- v1, v2: '(%d+),(%d+)[; ]+'
    for v1 in string.gfind(str, '-?(%d+)') do
        v1 = tonumber( v1 )
        if v1 < 1 or v1 > segCt then return false end
        if structure.GetAminoAcid(v1) ~= 'c' then return false end
        cysList[#cysList+1] = v1
    end
    if #cysList ~= 2*ctPairsWanted then return false end
    for i = 1, #cysList-1 do
        val = cysList[i]
        for j=i+1, #cysList do
            if cysList[j+1] == val then return false end
        end
    end

    -- ok, we believe the list is a correct number of distinct cysteines
    cysteinePairs = {}
    for i=1, ctPairsWanted do
        cysteinePairs[#cysteinePairs+1] = {cysList[2*i-1], cysList[2*i]}
    end
end

function SetCysteineOptions( dlg )
    ctCysPairs = dlg.cyPairsWanted.value
    if dlg.cyPairings.value ~= nil then
        parseCysteineStr( ctCysPairs, dlg.cyPairings.value)
    else
        -- this happens if user puts empty comment, 
        -- but says they want a number of pairs smaller than expected
        --  we interpret this as request to truncate accepted list
        while ctCysPairs < #cysteinePairs do
            table.remove(cysteinePairs, #cysteinePairs)
        end
    end
	
	if ctCysPairs > 0 then
		str = ""
		for i=1, #cysteinePairs do
			pair = cysteinePairs[i]
			str = str..pair[1]..","..pair[2]
			if i < #cysteinePairs then
				str = str .. "; "
			end
		end
		DebugPrint( "CysteinePairs: "..str )
		print( "Try to add "..ctCysPairs.." cysteine bands" )
	else
		print( "No cysteine bands wanted, none made" )
	end
end

function LetUserChooseParams( )
    dlg = dialog.CreateDialog( "Simulated Annealing And More" )
    dlg.cmt1 = dialog.AddLabel( "Pick one of following (topmost checked will be run)" )
    dlg.wantSASeq = dialog.AddCheckbox( "SA - Sequence", true )
    dlg.cmt2 = dialog.AddLabel( "--  Single Options" )
    dlg.wantSAEarly = dialog.AddCheckbox( "SA - Early", false )
    dlg.wantSaRebuilder = dialog.AddCheckbox( "SA - Rebuild", false )
    dlg.wantSAMid = dialog.AddCheckbox( "SA - Midgame", false )
    dlg.cmt3 = dialog.AddLabel( "-- Oddball Single Options" )
    dlg.wantSALowCI = dialog.AddCheckbox( "SA - LowCI", false )
    dlg.wantSALate = dialog.AddCheckbox( "SA - Late Game", false )
    if PuzzleHasContactMap then
        dlg.wantSAContact = dialog.AddCheckbox( "SA - ContactMap", false )
    end
    dlg.modLabel = dialog.AddLabel( "------------------Option(s)-----------------------" )
    if #cysteineSegs >= 2 then
        dlg.playCysteineGames = dialog.AddCheckbox( "Bond cysteine pairs", false )
    end
    if PuzzleHasContactMap then
        dlg.cmForceContactBettering = dialog.AddCheckbox( "Contact map tries to improve", CONTACTMAP_TAKEDROPS )
        dlg.cmInfo = dialog.AddLabel( "If tries to improve: Max score drop allowed" )
        dlg.cmMaxDrop = dialog.AddSlider( "  Drop", CONTACTMAP_MAX_DROP, 5, 100, 0 )
    end
    if PuzzleAllowsMutate then
        dlg.cmMutateOrShake = dialog.AddCheckbox( "Mutate instead of shake", true )
        dlg.cmAllowPointMutate = dialog.AddCheckbox( "Perform single-aa mutations", ALLOW_SPOT_MUTATION )
    end
    if PuzzleUsesFilters then
        dlg.cmKillFilters = dialog.AddCheckbox( "DisableFilters", KILLING_FILTERS )
    end
    dlg.preserveCuts = dialog.AddCheckbox( "Preserve user cuts", PRESERVE_USER_CUTS )
    if InitialBandCount > 0 then
        dlg.preserveBands = dialog.AddCheckbox( "Preserve user bands", true )
    end
    dlg.allowExpanding = dialog.AddCheckbox( "Allow Expansion", true )
    dlg.moreOptions = dialog.AddCheckbox( "Advanced options", false )
    dlg.ok = dialog.AddButton( "OK", 1 )
    dlg.cancel = dialog.AddButton( "Cancel", 0 )
    if dialog.Show( dlg ) == 0 then return false end

    -- PICK the Operation
    if dlg.wantSASeq.value then OperationChosen = "SASequence"
    elseif dlg.wantSALowCI.value then OperationChosen = "SALowCI"
    elseif dlg.wantSAEarly.value then OperationChosen = "SAEarly"
    elseif dlg.wantSaRebuilder.value then OperationChosen = "SARebuild"
    elseif dlg.wantSAMid.value then OperationChosen = "SAMidgame"
    elseif dlg.wantSALate.value then OperationChosen = "SALate"
    elseif PuzzleHasContactMap and dlg.wantSAContact.value then OperationChosen = "SAContact"
    else return false -- something went wrong...
    end

    -- PICK other info from initial dialog box
    if PuzzleHasContactMap then
        CONTACTMAP_TAKEDROPS = dlg.cmForceContactBettering.value
        if CONTACTMAP_TAKEDROPS then 
            print( "Score drops of ".. CONTACTMAP_MAX_DROP .. " accepted if ContactMap improves" ) 
        else 
            print( "No score drops permitted" ) 
        end
    end
    if PuzzleAllowsMutate then
        MUTATE_NOT_SHAKE = dlg.cmMutateOrShake.value
        ALLOW_SPOT_MUTATION = dlg.cmAllowPointMutate.value
    end
    if PuzzleUsesFilters then
        KILLING_FILTERS = dlg.cmKillFilters.value
    end

    if InitialBandCount > 0 then
        BANDS_PRESERVE_USER = dlg.preserveBands.value
        if not BANDS_PRESERVE_USER then
            InitialBandCount = 0
            DelBands( )
        end
        DebugPrint( "Preserve user bands? ".. BoolStr( PRESERVE_USER_BANDS ) )
    end

    PRESERVE_USER_CUTS = dlg.preserveCuts.value
    if PRESERVE_USER_CUTS then
        IDEALIZE_USE_CUTS = false
        REBUILD_USECUTS = false
        PROB_FUSEWIGGLE_CUTFUSE = 0.00 -- don't use wigglers that put cuts
        PROB_FUSEWIGGLE_LOCAL = 0.00  -- this one's too useless
        PROB_FUSEWIGGLE_BLUEFUSE = 1.0
    end
    DebugPrint( "Preserve user cuts? ".. BoolStr( PRESERVE_USER_CUTS ) )
        
    ALLOW_FUSE = false
    if OperationChosen == "SALate" then 
        ALLOW_FUSE = true
        DebugPrint( "Final Fuse Run? true")
    end
    
    ALLOW_EXPAND_FOR_ACTION = dlg.allowExpanding.value
    DebugPrint( "Expand? ".. BoolStr( ALLOW_EXPAND_FOR_ACTION))
    
    -- DO WE PLAY WITH CYSTEINES?
    if dlg.playCysteineGames.value == true then
        PERFORM_CYSTEINE_BANDING = true
        dlg2 = dialog.CreateDialog("Cysteine Pairings")
        PutCysteineOptionsInBox( dlg2 )
        if dialog.Show( dlg2 ) ~= 0 then
			if InitialBandCount > 0 then 
				print("Cannot honor user bands if cysteine bands are added")
				InitialBandCount = 0
				DelBands( )
			end
            SetCysteineOptions( dlg2 )
			PutCysteinePairBands( )
        end
    end

    -- ARE WE DONE?
    if not dlg.moreOptions.value then return true end

    -- NO: USER WANTS MORE OPTIONS, offer as relevant to OperationChosen 
    if OperationChosen == "SASequence" then
        -- also allow user to turn on fuses here
        dlg = dialog.CreateDialog( "SA Sequence Parameters" )
        PutSASequenceOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, true )
        if dialog.Show( dlg ) == 0 then return false end
        SetSASequenceOptions( dlg )
        SetBasicSpecialOptions( dlg, true, true )
    elseif OperationChosen == "SAEarly" then
        dlg = dialog.CreateDialog( "Early SA Options" )
        PutEarlyOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, false )
        if dialog.Show( dlg ) == 0 then return false end
        SetEarlyOptions( dlg )
        SetBasicSpecialOptions( dlg, true, false )
    elseif OperationChosen == "SARebuild" then
        dlg = dialog.CreateDialog( "Rebuild-SA Parameters" )
        PutRebuildOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, false, false )
        if dialog.Show( dlg ) == 0 then return false end
        SetRebuildOptions( dlg )
        SetBasicSpecialOptions( dlg, false, false )
    elseif OperationChosen == "SAMidgame" then
        dlg = dialog.CreateDialog( "Mid SA Options" )
        PutMidOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, false )
        if dialog.Show( dlg ) == 0 then return false end
        SetMidOptions( dlg )
        SetBasicSpecialOptions( dlg, true, false )
        DebugPrint( "probTarg="..PROB_TARGETED )
    elseif OperationChosen == "SALate" then
        ALLOW_FUSE = true
        dlg = dialog.CreateDialog( "Late SA Options" )
        PutLateOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, true )
        if dialog.Show( dlg ) == 0 then return false end
        SetLateOptions( dlg )
        SetBasicSpecialOptions( dlg, true, true )
        
    elseif OperationChosen == "SAContact" then
        -- not in normal sequence, even for contact puzzles
        dlg = dialog.CreateDialog( "ContactMap SA Options" )
        PutContactMapOptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, false )
        if dialog.Show( dlg ) == 0 then return false end
        SetContactMapOptions( dlg )
        SetBasicSpecialOptions( dlg, true, false )
    elseif OperationChosen == "SALowCI" then
        -- not in normal sequence
        dlg = dialog.CreateDialog( "LowCI SA Options" )
        PutLowCI20OptionsInBox( dlg )
        PutLowCI50OptionsInBox( dlg )
        AddBasicSpecialOptions( dlg, true, false )
        if dialog.Show( dlg ) == 0 then return false end
        SetLowCI20Options( dlg )
        SetLowCI50Options( dlg )
        SetBasicSpecialOptions( dlg, true, false )
    end
 
    return true
end

--------------------------------------------------------------------------
-- SETUP and CLEANUP
--------------------------------------------------------------------------
function CleanPuzzleState( )
    ExactSetCI( InitialClashImportance )
    selection.DeselectAll( )
    ResetFrozenness( )
    ResetCuts( )
    DelBands( true ) 
    ResetSecondaryStructures( 1, segCt )
end

function PrintState( )
    local gain = getScore( ) - InitialScore
    if gain < 0.001 then 
        print( " No change" )
    else
        print( " Startscore: "..TrimNum( InitialScore ))
        print( " Score: "..TrimNum( getScore( ) ))
        print( " Total gain: "..TrimNum( gain ) )
    end
    print( " Run time: "..os.time( ) - StartTime )
end

function End( errstr )
    if EndCalled then return end    -- no infinite recursion please
    EndCalled = true

    print( "" )
    if errstr ~= nil then
        if SHOW_STATS_REPORT then
            print("")
            ReportStats( )
            print( "TotalShakeTime="..TotalShakeTime )   -- because I want to know.
        end
        
        if string.find( errstr, "Cancelled" ) then
            print( "User cancel" )
        else
            print( errstr )
        end
    end
    PrintState( )

    if KILLING_FILTERS then ResetSlowFilterState( ) end
    if RESTORE_BEST_AT_END then LoadBest( ) end
    if PERFORM_CLEANUP_AT_END then CleanPuzzleState( ) end
end

function InitializePuzzleState( )
    seedRandom( )
    InitialScore = getScore( )
    save.Quicksave( QS_Start )
    save.Quicksave( QS_CmBest )
    save.Quicksave( QS_Best )
    StartTime = os.time( )

    -- learn a few things
    InitialClashImportance = behavior.GetClashImportance( )
    SCALE_MAXCI = InitialClashImportance
    PuzzleAllowsMutate = IsPuzzleMutable( )
    if PuzzleAllowsMutate then PuzzleIsDesigner = checkAla( ) end
    PuzzleUsesFilters = CheckForSlowFilters( )
    KILLING_FILTERS = false  -- PuzzleUsesFilters
    CheckForContactMap( )
    if PuzzleHasContactMap then InitializeContactMapSegList( CONTACTMAP_THRESHOLD ) end

    GenerateNonLockedSegList( )
    StoreSecondaryStructure( )
    SaveFrozenness( )
    FindAllMovableSegs( )
    LocateHelices( )
    LocateSheets( )
    LocateLoops( )
    LocateMissingHelices( )  -- places where a possible helix may be hiding
    if #possibleHelices > 0 then
        print( "Treating some loops as helices: count of ranges ="..#possibleHelices )
    end
    SetCysteineSegList( )
    GuessCysteinePairs( )
    DelBands( )
    selection.DeselectAll( )
    recentbest.Save( )  -- do this after DelBands
end

--------------------------------------------------------------------------
-- MAIN
--------------------------------------------------------------------------
function main( )
    OperationChosen = "none"
    InitializePuzzleState()
    didSomething = false
    print( "Start is in qs="..QS_Start.."; Best is in "..QS_Best )
    if PuzzleHasContactMap then
        print( "CM-best is in "..QS_CmBest )
    end
    
    ------------------------DEBUG STUFF
    -- if any actions could stand some debugging, put them here
    -- then multi-line comment from below to end of function
    ------------------------
    
    while LetUserChooseParams( ) do
        if KILLING_FILTERS then TurnOffSlowFilters( ) end
        PerformOperation( OperationChosen )
        if OperationChosen == "SASequence" then
            print("End of Early - qs 4")
            print("End of Rebuild - qs 5")
            if ALLOW_FUSE then
                print("End of Midgame - qs 6")
            end
        end
        didSomething = true
        LoadBest( )
        PrintState( )
        if SHOW_STATS_REPORT then
            ReportStats()
        end
        ResetStats()
    end
    
    print( " Done!" )
    if didSomething then End( ) end
end

xpcall( main, End )