tests/example.py

   1 import sys
   2 import pysam
   3
   4 samfile = pysam.Samfile( "ex1.bam", "rb" )
   5
   6 print "###################"
   7 # check different ways to iterate
   8 print len(list(samfile.fetch()))
   9 print len(list(samfile.fetch( "chr1", 10, 200 )))
  10 print len(list(samfile.fetch( region="chr1:10-200" )))
  11 print len(list(samfile.fetch( "chr1" )))
  12 print len(list(samfile.fetch( region="chr1")))
  13 print len(list(samfile.fetch( "chr2" )))
  14 print len(list(samfile.fetch( region="chr2")))
  15 print len(list(samfile.fetch()))
  16 print len(list(samfile.fetch( "chr1" )))
  17 print len(list(samfile.fetch( region="chr1")))
  18 print len(list(samfile.fetch()))
  19
  20 print len(list(samfile.pileup( "chr1", 10, 200 )))
  21 print len(list(samfile.pileup( region="chr1:10-200" )))
  22 print len(list(samfile.pileup( "chr1" )))
  23 print len(list(samfile.pileup( region="chr1")))
  24 print len(list(samfile.pileup( "chr2" )))
  25 print len(list(samfile.pileup( region="chr2")))
  26 print len(list(samfile.pileup()))
  27 print len(list(samfile.pileup()))
  28
  29 print "########### fetch with callback ################"
  30 def my_fetch_callback( alignment ): print str(alignment)
  31 samfile.fetch( region="chr1:10-200", callback=my_fetch_callback )
  32
  33 print "########## pileup with callback ################"
  34 def my_pileup_callback( column ): print str(column)
  35 samfile.pileup( region="chr1:10-200", callback=my_pileup_callback )
  36
  37 print "##########iterator row #################"
  38 iter = pysam.IteratorRow( samfile, 0, 10, 200)
  39 for x in iter: print str(x)
  40
  41 print "##########iterator col #################"
  42 iter = pysam.IteratorColumn( samfile, 0, 10, 200 )
  43 for x in iter: print str(x)
  44
  45 print "#########row all##################"
  46 iter = pysam.IteratorRowAll( samfile )
  47 for x in iter: print str(x)
  48
  49
  50 print "###################"
  51
  52 class Counter:
  53     mCounts = 0
  54     def __call__(self, alignment):
  55         self.mCounts += 1
  56
  57 c = Counter()
  58 samfile.fetch( "chr1:10-200", c )
  59 print "counts=", c.mCounts
  60
  61 sys.exit(0)
  62 print samfile.getTarget( 0 )
  63 print samfile.getTarget( 1 )
  64
  65 for p in pysam.pileup( "-c", "ex1.bam" ):
  66     print str(p)
  67
  68 print pysam.pileup.getMessages()
  69
  70 for p in pysam.pileup( "-c", "ex1.bam", raw=True ):
  71     print str(p),
  72
  73
  74
  75 print "###########################"
  76
  77 samfile = pysam.Samfile( "ex2.sam.gz", "r" )
  78
  79 print "num targets=", samfile.getNumTargets()
  80
  81 iter = pysam.IteratorRowAll( samfile )
  82 for x in iter: print str(x)
  83
  84 samfile.close()
  85
  86 print "###########################"
  87 samfile = pysam.Samfile( "ex2.sam.gz", "r" )
  88 def my_fetch_callback( alignment ):
  89     print str(alignment)
  90
  91 try:
  92     samfile.fetch( "chr1:10-20", my_fetch_callback )
  93 except AssertionError:
  94     print "caught fetch exception"
  95
  96 samfile.close()
  97
  98 print "###########################"
  99 samfile = pysam.Samfile( "ex2.sam.gz", "r" )
 100 def my_pileup_callback( pileups ):
 101     print str(pileups)
 102 try:
 103     samfile.pileup( "chr1:10-20", my_pileup_callback )
 104 except NotImplementedError:
 105     print "caught pileup exception"
 106
 107 # playing arount with headers
 108 samfile = pysam.Samfile( "ex3.sam", "r" )
 109 print samfile.targets
 110 print samfile.lengths
 111 print samfile.text
 112 print samdile.header
 113 header = samfile.header
 114 samfile.close()
 115
 116 header["HD"]["SO"] = "unsorted"
 117 outfile = pysam.Samfile( "out.sam", "wh",
 118                          header = header )
 119
 120 outfile.close()
 121