import sys
import teem
import pullDemo

# this returns a ctypes-allocated array of meetPullVols, which are a way of 
# storing volumes prior to being passed to pull.
# Format is <vol>:<kind>[:<minScale>-<#SclSamples>-<maxScale>[-o|n]]:<volName>
# where the flag "-o" means use pre-computed sigma locations optimized
# for accurate reconstruction, and "-n" means do scale-normalized 
# derivatives, which are necessary for meaningful cross-scale 
# feature strength comparisons
vol = pullDemo.volLoad(['bntorus-c4h.nrrd:scalar:0-11-10-o:V',
                        'bntorus-c4h.nrrd:scalar:0-11-10-on:VSN'],
                       'blur', 'ds:1,5', 1)

# How the information for the particle system is learned in the volumes
# via gage. Format is <info>:<volName>:<gageItem>[:zero:sign]
infoList = ['spthr:V:val:0.3:1',
            'sthr:VSN:heval1:-4:-1',
            'strn:VSN:heval1:0:-1',
            'h-c:V:val:0:-1',
            'hgvec:V:gvec',
            'hhess:V:hess',
            'tan1:V:hevec2',
            'tan2:V:hevec1']

Phi1 = {'type':teem.pullInterTypeUnivariate, 'r':'cwell:0.6,-0.002'}
Phi2 = {'type':teem.pullInterTypeAdditive, 'r':'cwell:0.6,-0.002',
        's':'bparab:20,0.88,-0.03', 'w':'butter:20,0.88'}

# create container for particle positions
npos = teem.nrrdNew()

radSpace = 0.07
lti = True
cbst = True
# Run system with Phi_1 energy: sampling full scale-space extent of feature
pullDemo.run(npos,
             init=[teem.pullInitMethodPointPerVoxel, -3, 1, 0, 4, 1],
             vol=vol,
             info=infoList + ['lthr:VSN:heval1:-5:-1'],
             energy=Phi1,
             cbst=cbst, efs=False, lti=lti, npcwza=False,
             radSpace=radSpace,
             radScale=1.2,
             alpha=0.5,
             beta=0,
             iterMax=100)

if (teem.nrrdSave("torus-phase1.nrrd", npos, None)):
    estr = teem.biffGetDone('nrrd')
    print "problem saving output:\n%s" % estr
    sys.exit(1)

# Run system with Phi_2 energy but with alpha=0, to congregate the 
# samples along the scale of maximal ridge strength

pullDemo.run(npos,
             init=[teem.pullInitMethodGivenPos, npos],
             vol=vol,
             info=infoList + ['lthr:VSN:heval1:-5:-1'],
             energy=Phi2,
             cbst=cbst, efs=True, lti=lti, npcwza=True,
             radSpace=radSpace,
             radScale=4,
             alpha=0.0,
             beta=0.25,
             gamma=1,
             iterMax=20)

if (teem.nrrdSave("torus-phase2.nrrd", npos, None)):
    estr = teem.biffGetDone('nrrd')
    print "problem saving output:\n%s" % estr
    sys.exit(1)
    
# Recover the gamma that was cleverly learned (by balling pullGammaLearn)
# at the end of the last run
gam = pullDemo.gamma
print 'Learned gamma = %g\n' % gam

# Final run to stabilize population along ridge. Uses a rather
# aggressive (200) live threshold ("lthr") on feature strength 
# so as to remove stragglers
pullDemo.run(npos,
             init=[teem.pullInitMethodGivenPos, npos],
             vol=vol,
             info=infoList + ['lthr:VSN:heval1:-200:-1'],
             energy=Phi2,
             cbst=cbst, efs=True, lti=lti, npcwza=True,
             radSpace=radSpace,
             radScale=4,
             alpha=0.25,
             beta=0.25,
             gamma=gam,
             edmin=0.0000001,
             iterMax=220)

# Free the loaded volumes
pullDemo.volFree(vol)

# Save results
if (teem.nrrdSave("torus-phase3.nrrd", npos, None)):
    estr = teem.biffGetDone('nrrd')
    print "problem saving output:\n%s" % estr
    sys.exit(1)