Python Example Scripts

See Python Interface for introduction to running phaser as a python library.

Example scripts for the most popular modes of running Phaser.

Reading MTZ Files for Molecular Replacement

Example script for reading data from MTZ file beta_blip.mtz.

Note that by default reflections are sorted into resolution order upon reading, to achieve a performance gain in the molecular replacement routines. If reflections are not being read from an MTZ file with this script, reflections should be pre-sorted into resolution order to achieve the same performance gain. Sorting is turned off with the setSORT(False) function.

     #beta_blip_data.py
     from phaser import *
     i = InputMR_DAT()
     HKLIN = "beta_blip_P3221.mtz"
     i.setHKLI(HKLIN)
     i.setMUTE(True)
     r = runMR_DAT(i)
     print r.logfile()
     if r.Success():
       hkl = r.getMiller()
       fobs = r.getF()
       sigma = r.getSIGF()
       nrefl = min(10,hkl.size())
       print "Data read from: " , HKLIN
       print "First ", nrefl , " reflections"
       print "%4s %4s %4s %10s %10s" % ("H","K","L",F,SIGF)
       for i in range(0,nrefl):
         print "%4d %4d %4d %10.4f %10.4f" % \
           (hkl[i][0],hkl[i][1],hkl[i][2],fobs[i],sigma[i])
     else:
       print "Job exit status FAILURE"
       print r.ErrorName(), "ERROR :", r.ErrorMessage()

Automated Molecular Replacement

Example script for automated structure solution of BETA-BLIP

     #beta_blip_auto.py
     from phaser import *
     i = InputMR_DAT()
     i.setHKLI("beta_blip_P3221.mtz")
     i.setHIRES(6.0)
     i.setMUTE(True)
     r = runMR_DAT(i)
     if r.Success():
       i = InputMR_AUTO()
       i.setREFL_DATA(r.getREFL_DATA())
       i.setROOT("beta_blip_auto")
       i.addENSE_PDB_ID("beta","beta.pdb",1.0)
       i.addENSE_PDB_ID("blip","blip.pdb",1.0)
       i.addCOMP_PROT_MW_NUM(28853,1)
       i.addCOMP_PROT_MW_NUM(17522,1)
       i.addSEAR_ENSE_NUM("beta",1)
       i.addSEAR_ENSE_NUM("blip",1)
       i.setMUTE(True)
       del(r)
       r = runMR_AUTO(i)
       if r.Success():
         if r.foundSolutions() :
           print "Phaser has found MR solutions"
           print "Top LLG = %f" % r.getTopLLG()
           print "Top PDB file = %s" % r.getTopPdbFile()
         else:
           print "Phaser has not found any MR solutions"
       else:
         print "Job exit status FAILURE"
         print r.ErrorName(), "ERROR :", r.ErrorMessage()
     else:
       print "Job exit status FAILURE"
       print r.ErrorName(), "ERROR :", r.ErrorMessage()

Reading MTZ Files for Experimental Phasing

Example script for reading SAD data from MTZ file S-insulin.mtz

      #insulin_data.py
      from phaser import *
      i = InputEP_DAT()
      HKLIN = "S-insulin.mtz"
      xtalid = "insulin"
      waveid = "cuka"
      i.setHKLI(HKLIN)
      i.addCRYS_ANOM_LABI(xtalid,waveid,"F(+)","SIGF(+)","F(-)","SIGF(-)")
      i.setMUTE(True)
      r = runEP_DAT(i)
      if r.Success():
        hkl = r.getMiller()
        Fpos = r.getFpos(xtalid,waveid)
        Ppos = r.getPpos(xtalid,waveid)
        Fneg = r.getFneg(xtalid,waveid)
        Pneg = r.getPneg(xtalid,waveid)
        print "Data read from: " , HKLIN
        print "Spacegroup Name (Hall symbol) = %s (%s)" % \
         (r.getSpaceGroupName(), r.getSpaceGroupHall())
        print "Unitcell = " , r.getUnitCell()
        nrefl = min(10,hkl.size())
        print "First ", nrefl , " reflections with anomalous differences"
        print "%4s %4s %4s %10s %10s %10s" % ("H","K","L","F(+)","F(-)","D")
        i = 0
        r = 0
        while r < nrefl:
          if Ppos[i] and Pneg[i] :
            D = abs(Fpos[i]-Fneg[i])
            if D > 0
              print "%4d %4d %4d %10.4f %10.4f %10.4f" % \
                (hkl[i][0],hkl[i][1],hkl[i][2],Fpos[i],Fneg[i],D)
              r=r+1
          i=i+1
      else:
        print "Job exit status FAILURE"
        print r.ErrorName(), "ERROR :", r.ErrorMessage()

Automated Experimental Phasing

Example script for SAD phasing for insulin

      #insulin_sad.py
      from phaser import *
      from cctbx import xray
      i = InputEP_DAT()
      HKLIN = "S-insulin.mtz"
      xtalid = "insulin"
      waveid = "cuka"
      i.setHKLI(HKLIN)
      i.addCRYS_ANOM_LABI(xtalid,waveid,"F(+)","SIGF(+)","F(-)","SIGF(-)")
      i.setMUTE(True)
      r = runEP_DAT(i)
      if r.Success():
        hkl = r.getMiller()
        Fpos = r.getFpos(xtalid,waveid)
        Spos = r.getSIGFpos(xtalid,waveid)
        Ppos = r.getPpos(xtalid,waveid)
        Fneg = r.getFneg(xtalid,waveid)
        Sneg = r.getSIGFneg(xtalid,waveid)
        Pneg = r.getPneg(xtalid,waveid)
        i = InputEP_AUTO()
        i.setSPAC_HALL(r.getSpaceGroupHall())
        i.setCELL(r.getUnitCell())
        i.setCRYS_MILLER(hkl)
        i.addCRYS_ANOM_DATA(xtalid,waveid,Fpos,Spos,Ppos,Fneg,Sneg,Pneg)
        i.setATOM_PDB(xtalid,"S-insulin_hyss.pdb")
        i.setLLGC_CRYS_COMPLETE(xtalid,True)
        i.addLLGC_CRYS_SCAT_ELEMENT(xtalid,"S")
        i.addCOMP_PROT_FASTA_NUM("S-insulin.seq",1.)
        i.setHKLO(False)
        i.setSCRI(False)
        i.setXYZO(False)
        i.setMUTE(True)
        r = runEP_AUTO(i)
        if r.Success():
          print "SAD phasing"
          print "Data read from: " , HKLIN
          print "Data output to : " , r.getMtzFile()
          print "Spacegroup Name (Hall symbol) = %s (%s)" % \
            (r.getSpaceGroupName(), r.getSpaceGroupHall())
          print "Unitcell = " , r.getUnitCell()
          print "LogLikelihood = " , r.getLogLikelihood()
          atom = r.getAtoms(xtalid)
          print atom.size(), " refined atoms"
          print "%5s %10s %10s %10s %10s %10s" % \
            ("atom","x","y","z","occupancy","u-iso")
          for i in range(0,atom.size()):
            print "%5s %10.4f %10.4f %10.4f %10.4f %10.4f" % \
           
(atom[i].scattering_type,atom[i].site[0],atom[i].site[1],atom[i].site[2],atom[i].occupancy,atom[i].u_iso)
          hkl = r.getMiller();
          fwt = r.getFWT()
          phwt = r.getPHWT()
          fom = r.getFOM()
          nrefl = min(10,hkl.size())
          print "First ", nrefl , " reflections"
          print "%4s %4s %4s %10s %10s %10s" % \
            ("H","K","L","FWT","PHWT","FOM")
          for i in range(0,nrefl):
            print "%4d %4d %4d %10.4f %10.4f %10.4f" % \
              (hkl[i][0],hkl[i][1],hkl[i][2],fwt[i],phwt[i],fom[i])
        else:
          print "Job exit status FAILURE"
          print r.ErrorName(), "ERROR :", r.ErrorMessage()
      else:
        print "Job exit status FAILURE"
        print r.ErrorName(), "ERROR :", r.ErrorMessage()

Anisotropy Correction

Example script script for anisotropy correction of BETA-BLIP data

     #beta_blip_ano.py
     from phaser import *
     i = InputMR_DAT()
     HKLIn = "beta_blip_P3221.mtz"
     i.setHKLI(HKLIn)
     i.setMUTE(True)
     r = runMR_DAT(i)
     if r.Success():
       i = InputANO()
       i.setSPAC_HALL(r.getSpaceGroupHall())
       i.setREFL_DATA(r.getREFL_DATA())
       i.setROOT("beta_blip_ano")
       i.setMUTE(True)
       del(r)
       r = runANO(i)
       if r.Success():
         print "Anisotropy Correction"
         print "Data read from: " , HKLIn
         print "Data output to : " , r.getMtzFile()
         print "Spacegroup Name (Hall symbol) = %s (%s)" % \
           (r.getSpaceGroupName(), r.getSpaceGroupHall())
         print "Unitcell = " , r.getUnitCell()
         print "Principal components = " , r.getEigenBs()
         print "Range of principal components = " , r.getAnisoDeltaB()
         print "Wilson Scale = " , r.getWilsonK()
         print "Wilson B-factor = " , r.getWilsonB()
         hkl = r.getMiller();
         f = r.getF()
         sigf = r.getSIGF()
         f_iso = r.getCorrectedF(True)
         sigf_iso = r.getCorrectedSIGF(True)
         corr = r.getCorrection(True)
         nrefl = min(10,hkl.size())
         print "First ", nrefl , " reflections"
         print "%4s %4s %4s %10s %10s %10s %10s %10s" % \
           ("H","K","L","F","SIGF",r.getLaboutF(),r.getLaboutSIGF(),"Corr\'n")
         for i in range(0,nrefl):
           print "%4d %4d %4d %10.4f %10.4f %10.4f %10.4f %10.4f" % \
             (hkl[i][0],hkl[i][1],hkl[i][2],f[i],sigf[i],f_iso[i],sigf_iso[i],corr[i])
       else:
         print "Job exit status FAILURE"
         print r.ErrorName(), "ERROR :", r.ErrorMessage()
     else:
       print "Job exit status FAILURE"
       print r.ErrorName(), "ERROR :", r.ErrorMessage()

Cell Content Analysis

Example script for cell content analysis of BETA-BLIP

      #beta_blip_cca.py
      from phaser import *
      i = InputMR_DAT()
      HKLIN = "beta_blip.mtz"
      F = "Fobs"
      SIGF = "Sigma"
      i.setHKLI(HKLIN)
      i.setLABI_F_SIGF(F,SIGF)
      i.setMUTE(True)
      r = runMR_DAT(i)
      if r.Success():
        i = InputCCA()
        i.setSPAC_HALL(r.getSpaceGroupHall())
        i.setCELL6(r.getUnitCell())
        i.addCOMP_PROT_MW_NUM(28853,1)
        i.addCOMP_PROT_MW_NUM(17522,1)
        i.setMUTE(True)
        del(r)
        r = runCCA(i)
        if r.Success():
          print "Cell Content Analysis"
          print "Molecular weight of assembly = " , r.getAssemblyMW()
          print "Best Z value = " , r.getBestZ()
          print "Best VM value = " , r.getBestVM()
          print "Probability of Best VM = " , r.getBestProb()
        else:
          print "Job exit status FAILURE"
          print r.ErrorName(), "ERROR :", r.ErrorMessage()
      else:
        print "Job exit status FAILURE"
        print r.ErrorName(), "ERROR :", r.ErrorMessage()

Translational NCS Analysis

Example script for translational NCS analysis of BETA-BLIP

      #beta_blip_ncs.py
      from phaser import *
      i = InputMR_DAT()
      HKLIN = "beta_blip.mtz"
      F = "Fobs"
      SIGF = "Sigma"
      i.setHKLI(HKLIN)
      i.setLABI_F_SIGF(F,SIGF)
      i.setMUTE(True)
      r = runMR_DAT(i)
      if r.Success():
        i = InputNCS()
        i.setSPAC_HALL(r.getSpaceGroupHall())
        i.setCELL6(r.getUnitCell())
        i.setREFL_F_SIGF(r.getMiller(),r.getF(),r.getSIGF())
        i.addCOMP_PROT_MW_NUM(28853,1)
        i.addCOMP_PROT_MW_NUM(17522,1)
        i.setMUTE(True)
        del(r)
        r = runNCS(i)
        if r.Success():
          print "Translational NCS analysis"
          print "Translational NCS present = ", r.hasTNCS()
          if r.hasTNCS():
            print "Translational NCS vecor = ", r.hasTNCS()
          print "Twinning alpha = ", r.getTwinAlpha()
        else:
          print "Job exit status FAILURE"
          print r.ErrorName(), "ERROR :", r.ErrorMessage()
      else:
        print "Job exit status FAILURE"
        print r.ErrorName(), "ERROR :", r.ErrorMessage()

Normal Mode Analysis

Example script for normal mode analysis of BETA-BLIP. Note that the space group and unit cell are not required, and so the MTZ file does not need to be read to extract these parameters.

      #beta_nma.py
      from phaser import *
      i = InputNMA()
      i.setROOT("beta_nma")
      i.addENSE_PDB_ID("beta","beta.pdb",1.0)
      i.addNMA_MODE(4)
      i.setPERT_RMS_DIRE("FORWARD")
      i.setMUTE(True)
      r = runNMAXYZ(i)
      if r.Success():
        print "Normal Mode Analysis"
        for i in range(0,r.getNum()):
          print "PDB file = ", r.getPdbFile(i)
          displacement = r.getDisplacements(i)
          mode = r.getModes(i)
          for j in range(0,mode.size()):
            print " Mode = " , mode[j], " Displacement = ", displacement[j]
      else:
        print "Job exit status FAILURE"
        print r.ErrorName(), "ERROR :", r.ErrorMessage()

Logfile Handling

Example of how to redirect phaser output to a python string for real-time viewing of output, but not via standard output. Output to standard out is silenced with setMUTE(True).

      beta_blip_logfile.py
      from phaser import *
      from cStringIO import StringIO
      i = InputMR_DAT()
      i.setHKLI("beta_blip.mtz")
      i.setLABI_F_SIGF("Fobs","Sigma")
      i.setMUTE(True)
      o = Output()
      redirect_str = StringIO()
      o.setPackagePhenix(file_object=redirect_str)
      r = runMR_DAT(i,o)