def loadpair(file):
	id = file.split('/')[-1][:4]
	loadcenters(file,id+'nat')
	loadcenters(file[:-17]+'_decoy_sasa_centers.pdb',id+'decoy')
	d = id+'decoy'
	n = id+'nat'
	cmd.align(n,d)

    def testSaveRef(self, format):
        # for rms_cur (not all formats save all identifiers)
        m = -1
        cmd.set('retain_order')

        cmd.fragment('ala', 'm1')
        cmd.copy('m2', 'm1')
        cmd.copy('m3', 'm1')

        cmd.rotate('y', 90, 'm2')
        cmd.align('m3', 'm2')

        # with ref=m3
        with testing.mktemp('.' + format) as filename:
            cmd.save(filename, 'm2', ref='m3')
            cmd.load(filename, 'm4')
            self.assertAlmostEqual(cmd.rms_cur('m4', 'm1', matchmaker=m), 0.00, delta=1e-2)
            self.assertAlmostEqual(cmd.rms_cur('m4', 'm2', matchmaker=m), 1.87, delta=1e-2)

        # without ref
        with testing.mktemp('.' + format) as filename:
            cmd.save(filename, 'm2')
            cmd.load(filename, 'm5')
            self.assertAlmostEqual(cmd.rms_cur('m5', 'm2', matchmaker=m), 0.00, delta=1e-2)
            self.assertAlmostEqual(cmd.rms_cur('m5', 'm1', matchmaker=m), 1.87, delta=1e-2)

def sewalign(first,second):
	#hide waters
	cmd.hide("(solvent and (all))")
	#hide hydrogens
	cmd.hide( 'sticks', 'elem H' )
	cmd.hide( 'lines', 'elem H' )

	#show cartoon
	cmd.show("cartoon"   ,"all")

	#create duplicate of first
	first_copy = first + "_copy"
	cmd.copy(first_copy, first)

	#select first 14 residues 
	cmd.select('node_9_selection', '(obj *_9_* and resi 1-14)')
	cmd.select('node_12_selection', '(obj *_12_* and resi 1-14)')
	cmd.select('node_15_selection', '(obj *_15_* and resi 1-14)')

	alignment_1 = cmd.align(first, 'node_9_selection')
	print alignment_1[0]

	alignment_2 = cmd.align(second, 'node_12_selection')
	print alignment_2[0]

	alignment_3 = cmd.align(first_copy, 'node_15_selection')
	print alignment_3[0]

	def revert(self):
		v = cmd.get_view()
		cmd.remove(self.rdes.obj+" and not chain A")
		m = cmd.get_model(self.rdes.obj)
		n = self.oldcrd
		if not n:
			cmd.create("tmp12345",self.rnat.sel())
			cmd.align("tmp12345",self.rdes.sel())
			n = cmd.get_model("tmp12345").atom
			cmd.delete("tmp12345")
		di,ni = 0,0
		while m.atom[di].resi != str(self.rdes.resi): di += 1
		dj = di
		while m.atom[dj].resi == str(self.rdes.resi): dj += 1
		if self.oldcrd: self.oldcrd = None
		else:           self.oldcrd = m.atom[di:dj]
		m.atom = m.atom[:di] + n + m.atom[dj:]
		for i in range(di,di+len(n)):
			m.atom[i].resi  = str(self.rdes.resi)
			m.atom[i].chain = str(self.rdes.chain)
		cmd.load_model(m,self.rdes.obj,1)
		cmd.save("tmp.pdb",self.rdes.obj)
		cmd.delete(self.rdes.obj)
		cmd.load("tmp.pdb",self.rdes.obj,1)
		cmd.show('car',self.rdes.obj)
		cmd.show('lines')		
		redopent(self.rdes.obj)
		cmd.set_view(v)
		print "revert removing "+"".join(self.aas)+" from aas!"
		self.aas = [getaa('A',self.rdes.resi,self.manager.d.obj)]

def align_and_color(moving, fixed, maximum=None, native_color="blue"):
    if maximum == "0":
        maximum = None
    alignment_name = "alignment_" + moving
    cmd.align(moving, fixed, object=alignment_name, cycles=0)
    cmd.disable(alignment_name)
    rmsd_b(moving + " & " + alignment_name, fixed + " & " + alignment_name)
    cmd.color(native_color, fixed)
    cmd.spectrum("b", "green_red", moving, 0, maximum)

 def testSaveAln(self):
     cmd.fab('ACDEFGH', 'm1')
     cmd.fab('ACDFGH', 'm2')
     cmd.align('m1', 'm2', cycles=0, object='aln')
     with testing.mktemp('.aln') as filename:
         cmd.save(filename)
         with open(filename) as handle:
             lines = list(handle)
     self.assertEqual(lines[0].split(), ['CLUSTAL'])
     self.assertEqual(lines[1].split(), [])
     self.assertEqual(lines[2].split(), ['m1', 'ACDEFGH'])
     self.assertEqual(lines[3].split(), ['m2', 'ACD-FGH'])

def loadPackingPDB(file,name=None,native=None):
	"""
	usage: loadPackingPDB <file> , [<name for object>]
	loads a foo_packing.pdb file and colors it all pretty-like
	creates two selections along with the loaded object called
	NAMEcavities and NAMEprotein which are the heteratoms representing
	holes and everything else, respectively. Names can get pretty long,
	by pymol lets you do good stuff like "select NA*cav*", which will
	match a selection called NAMEISREALLYLONGcavities.
	"""

	if name is None:
		name = name = os.path.basename(file)
	if name.endswith('.gz'):
		name = name[:-3] 
	if name.endswith('.pdb'):
		name = name[:-4] 
	if name.endswith('.'):
		name = name[:-1] 
	if name.endswith("_packing"):
		name = name[:-8]

	zload(file,name)
	cmd.hide('everything',name)
	
	if native is not None:
		cmd.align(name,native)	
		cmd.zoom(native)

	useRosettaRadii()

	cavselname  = name+"cavities"
	protselname = name+"protein"

	cmd.select(cavselname, "resn CAV and b > 0.1 and %s"%(name) )
	cmd.select(protselname,"(not resn CAV) and %s"%(name) )

	useTempRadii(cavselname)
	useOccColors(cavselname)
	cmd.color("white",protselname)

	cmd.show('spheres', cavselname )
	cmd.show("cartoon",protselname)
	cmd.show("lines",protselname)

	cmd.select("none")
	cmd.delete("sele*")
	cmd.move('z',-50)

	return name

	def showpack(self):
		print self.toggle['sph']
		#cmd.hide('spheres','not '+self.manager.d.bupsel)
		cmd.hide('spheres')
		if self.toggle['sph']:
			self.toggle['sph'] = False
			return
		fr = self.rnat.obj+" and name n+ca+c and resi %i-%i"%(self.rnat.resi-10,self.rnat.resi+10)		
		to = self.rdes.obj+" and name n+ca+c and resi %i-%i"%(self.rdes.resi-10,self.rdes.resi+10)
		cmd.align(fr,to)
		cmd.center(self.rdes.sel()+" or "+self.rnat.sel())
		#cmd.show('sticks',self.rnat.sel())
		cmd.show('spheres',"byres ((not "+self.rnat.obj+") within 6 of "+self.rdes.sel()+")")
		self.toggle['sph'] = True

	def focus(self):
		self.manager.m = self
		if self.manager.prevm: cmd.hide('sticks',self.manager.prevm.rdes.sel())
		if self.manager.prevm: cmd.color('green',self.manager.prevm.rdes.sel()+" and elem C")
		fr = self.rnat.obj+" and name n+ca+c and resi %i-%i"%(self.rnat.resi-10,self.rnat.resi+10)		
		to = self.rdes.obj+" and name n+ca+c and resi %i-%i"%(self.rdes.resi-10,self.rdes.resi+10)
		cmd.align(fr,to)
		cmd.show('sticks',self.rdes.sel())
		cmd.color('white',self.rdes.sel()+" and elem C")
		cmd.center(self.rdes.sel())
		v = list(cmd.get_view())
		v[11] = -80.0
		cmd.set_view(tuple(v))
		cmd.clip('slab',50,self.rdes.sel())

    def testAlignMissingCoords(self):
        filename = self.datafile('1t46-frag.pdb')

        cmd.load(filename, 'm1')
        cmd.remove('m1 & resi 600-605')
        cmd.load(filename, 'm1', state=2)

        cmd.load(filename, 'm2')
        cmd.remove('m2 & resi 620-625')
        cmd.load(filename, 'm2', state=2)

        cmd.align('m1', 'm2', object='aln', mobile_state=1, target_state=1, cycles=0)

        self.assertIn('aln', cmd.get_names())

 def testSaveAlnNucleic(self):
     cmd.load(self.datafile('1rna.cif'))
     cmd.create('m1', 'chain A and not resi 6-7')
     cmd.create('m2', 'chain B')
     cmd.alter('m1 and resi 1', 'resn = "DT"') # mimic DNA
     cmd.alter('m2 and resi 20', 'resn = "UNK"') # mimic nonstd residue
     cmd.align('m1', 'm2', cycles=0, object='aln')
     with testing.mktemp('.aln') as filename:
         cmd.save(filename)
         with open(filename) as handle:
             lines = list(handle)
     self.assertEqual(lines[0].split(), ['CLUSTAL'])
     self.assertEqual(lines[1].split(), [])
     self.assertEqual(lines[2].split(), ['m1', 'TUAUA--UAUAUAA'])
     self.assertEqual(lines[3].split(), ['m2', 'UUAUA?AUAUAUAA'])

def ie_build_file(fname, align='true', ortho='true', hide='true',
                  zoom='true', col='[0.5, 0.5, 0.5]', scale='1'):
    if to_bool(ortho):
        cmd.set('orthoscopic', 'true')
    cmd.load(fname)
    if to_bool(align):
        object_list = cmd.get_names()
        target = object_list.pop()
        for obj in object_list:
            cmd.align(obj, target)
    if to_bool(hide):
        cmd.hide('everything')
    ie_build_all(col, scale)
    if to_bool(zoom):
        cmd.zoom()

def prepare_data(arg1):
    cmd.load("%s_r_b.pdb" % arg1, zoom=0)
    cmd.load("%s_l_b.pdb" % arg1, zoom=0)
    cmd.load("%s_r_u.pdb" % arg1, zoom=0)
    cmd.load("%s_l_u.pdb" % arg1, zoom=0)

    save_surface("%s_r_b" % arg1)
    save_surface("%s_l_b" % arg1)
    save_surface("%s_r_u" % arg1)
    save_surface("%s_l_u" % arg1)

    cmd.align("%s_r_b" % arg1, "%s_r_u" % arg1, object="ra")
    cmd.align("%s_l_b" % arg1, "%s_l_u" % arg1, object="la")
    cmd.save("receptor.aln", "ra")
    cmd.save("ligand.aln", "la")

def sa(intra=False,rainbow=True):
  """
  Superimpose all open models onto the first one.
  This may not work well with selections.
  Option intra can be set to True to enable intra_fit first, for working with multi-state (nmr) pdbs.
  [Thanks to Kyle Beauchamp for this one]
  """
  AllObj=cmd.get_names("all")
  for x in AllObj:
    print(AllObj[0],x)
    if intra==True:
      cmd.intra_fit(x)
    if rainbow==True:
      cmd.util.chainbow(x)
    cmd.align(x,AllObj[0])
    cmd.zoom()

def getnative():
    v = cmd.get_view()
    nats = []
    for obj in cmd.get_object_list():
        if len(obj) == 4:
            nats.append(obj)
            continue
        pid = obj[:4]
        if pid in nats:
            continue
        print "fetching native", pid
        cmd.fetch(pid)
        cmd.remove(pid + " and not chain A")
        cmd.remove(pid + " and resn HOH")
        cmd.align(pid, obj)
    cmd.set_view(v)

def alignallrms(sele):
    r = {}
    for i in cmd.get_object_list():
        r[i] = cmd.align("fr52re", i)[0]
    for k, v in r.items():
        if v < 2:
            print k, v

	def load(self):
		self.manager.d = self
		self.manager.m = None
		self.manager.m = None # to also set prevm = None		
		# cmd.delete(self.manager.prevd.obj)
		# cmd.delete(self.manager.prevd.pid)
		# print self.manager.d.obj
		# print self.manager.prevd.obj		
		# sys.exit()
		cmd.delete("all")
		cmd.load(self.getloadname(),self.obj,1)
		cmd.remove(self.obj+" and not chain A")
		cmd.fetch(self.pid)
		cmd.remove("het or hydro")
		cmd.color('gray',self.pid+' and elem C')
		self.scores['rms'] = cmd.align(self.pid,self.obj+" and chain A")[0]
		cmd.hide('ev')
		cmd.show('lines')
		cmd.show("car")
		#cmd.disable(self.pid)
		redopent(self.obj)
		self.recalc_design_pos()
		self.read_data_dir('avg_deg')
		self.read_data_dir('ddG')
		self.read_data_dir('rot_boltz')
		self.read_bup()
		self.read_tot_scores()
		self.read_res_scores()
		cmd.orient(self.obj+" or pnt*")
		self.manager.m = self.remembermm
		if self.manager.m: self.manager.m.focus()
		if self.remembermv: cmd.set_view(self.remembermv)

def align_all( subset = [] ):
  """
  Superimpose all open models onto the first one.
  This may not work well with selections.

  This function is probably taken from https://daslab.stanford.edu/site_data/docs_pymol_rhiju.pdf
  """
  print("""This returns a list with 7 items:

    RMSD after refinement
    Number of aligned atoms after refinement
    Number of refinement cycles
    RMSD before refinement
    Number of aligned atoms before refinement
    Raw alignment score
    Number of residues aligned """)

  AllObj=cmd.get_names("all")
  for x in AllObj[1:]:
    #print(AllObj[0],x)
    subset_tag = ''
    if isinstance( subset, int ):
      subset_tag = ' and resi %d' % subset
    elif isinstance( subset, list ) and len( subset ) > 0:
      subset_tag = ' and resi %d' % (subset[0])
      for m in range( 1,len(subset)): subset_tag += '+%d' % subset[m]
    elif isinstance( subset, str ) and len( subset ) > 0:
      subset_tag = ' and %s' % subset
    values = cmd.align(x+subset_tag,AllObj[0]+subset_tag)
    print(AllObj[0], x, ' '.join([str(v) for v in values]), '-- RMSD', values[3], ' of ', values[6], 'residues')
    cmd.zoom()

def align_all(target=None,mobile_selection='name ca',target_selection='name ca',cutoff=2, cycles=5,cgo_object=0):
  """
  Aligns all models in a list to one target

  usage:
    align_all [target][target_selection=name ca][mobile_selection=name ca][cutoff=2][cycles=5][cgo_object=0]
        where target specifies the model id you want to align all others against,
        and target_selection, mobile_selection, cutoff and cycles are options
        passed to the align command.

    By default the selection is all C-alpha atoms and the cutoff is 2 and the
    number of cycles is 5.
    Setting cgo_object to 1, will cause the generation of an alignment object for
    each object.  They will be named like <object>_on_<target>, where <object> and
    <target> will be replaced by the real object and target names.

    Example:
      align_all target=name1, mobile_selection=c. b & n. n+ca+c+o,target_selection=c. a & n. n+ca+c+o

  """
  cutoff = int(cutoff)
  cycles = int(cycles)
  cgo_object = int(cgo_object)

  object_list = cmd.get_names()
  object_list.remove(target)

  rmsd = {}
  rmsd_list = []
  for i in range(len(object_list)):
    if cgo_object:
      objectname = 'align_%s_on_%s' % (object_list[i],target)
      rms = cmd.align('%s & %s'%(object_list[i],mobile_selection),'%s & %s'%(target,target_selection),cutoff=cutoff,cycles=cycles,object=objectname)
    else:
      rms = cmd.align('%s & %s'%(object_list[i],mobile_selection),'%s & %s'%(target,target_selection),cutoff=cutoff,cycles=cycles)

    rmsd[object_list[i]] = rms[0]
    rmsd_list.append((rms[0],object_list[i]))

  rmsd_list.sort()
# loop over dictionary and print out matrix of final rms values
  print "Aligning against:",target
  for object_name in object_list:
    print "%s: %6.3f" % (object_name,rmsd[object_name])

  for r in rmsd_list:
    print "%6.3f  %s" % r

def redoA(sel="not sub", N=None):
    cmd.delete("sub*")
    cmd.hide("ev", sel + " and not chain A")
    v = cmd.get_view()
    chains = list(set(a.chain for a in cmd.get_model("name ca").atom))
    chains.sort()
    if N:
        chains = chains[:N]
    for c in chains:
        if c == "A":
            continue
        cmd.create("sub" + c, sel + " and chain A")
        cmd.align("sub" + c, sel + " and chain " + c)
        cmd.alter("sub" + c, "chain = '%s'" % c)
    color_by_chain()
    cmd.set_view(v)
    print charge(sel + " and chain A")