def build_bowtie_index(self):
     io.execute([
             'bowtie2-build',
             self.reference_fasta_filename(),
             self/'bowtie',
             ],
         )

 def build_snpeff(self):
     jar = io.find_jar('snpEff.jar')
     
     with open(self/'snpeff.config','wb') as f:
         print >> f, 'data_dir = snpeff'
         print >> f, 'genomes : ' + self.name
         print >> f, self.name + '.genome : ' + self.name 
     
     snpwork = io.Workspace(self/'snpeff',must_exist=False)
     snpwork_genome = io.Workspace(snpwork/self.name,must_exist=False)
     snpwork_genomes = io.Workspace(snpwork/'genomes',must_exist=False)
     
     annotations = self.annotations_filename()
     assert annotations
     with open(snpwork_genome/'genes.gff','wb') as f:
         for record in annotation.read_annotations(annotations):
             if record.end <= record.start: continue
             if not record.attr:
                 record.attr['attributes'] = 'none'
             print >> f, record.as_gff()
     
     with open(snpwork_genomes/(self.name+'.fa'),'wb') as f:
         for name, seq in io.read_sequences(self.reference_fasta_filename()):
             io.write_fasta(f, name, seq)
             
     io.execute('java -jar JAR build NAME -gff3 -c CONFIG',
         JAR=jar, NAME=self.name, CONFIG=self/'snpeff.config')

 def run(self):
     genome = self.genome
     if os.path.isdir(genome):
         genome = os.path.join(genome, os.path.split(genome)[1]+'.genome')
         print genome
     
     #pref_filename = os.path.join(os.path.expanduser('~'),'igv','prefs.properties')
     #if os.path.exists(pref_filename):
     #    with open(pref_filename,'rb') as f:
     #        lines = f.readlines()
     #    with open(pref_filename,'wb') as f:
     #        for line in lines:
     #            if line.startswith('DEFAULT_GENOME_KEY='):
     #                #line = 'DEFAULT_GENOME_KEY=\n'
     #                continue
     #            f.write(line)
     
     with workspace.tempspace() as temp:
         with open(temp/'batch.txt','wb') as f:
             print >> f, 'new'
             print >> f, 'preference LAST_TRACK_DIRECTORY', os.getcwd()
             print >> f, 'preference LAST_GENOME_IMPORT_DIRECTORY', os.getcwd()
             print >> f, 'genome '+os.path.abspath(genome)
             for filename in self.files:
                 print >> f, 'load '+os.path.abspath(filename)
         
         io.execute(['java','-Xmx32000m',
                     #Flags from igb.sh script:
                     '-Dproduction=true','-Dapple.laf.useScreenMenuBar=true','-Djava.net.preferIPv4Stack=true',
                     '-jar',io.find_jar('igv.jar'),'-b',temp/'batch.txt'])

def sort_and_index(in_filename, out_prefix):
    io.execute([
        'samtools', 'sort', in_filename, out_prefix
    ])
    
    io.execute([
        'samtools', 'index', out_prefix + '.bam'
    ])

 def build_bowtie_index(self, log_to=sys.stdout):
     io.execute([
             'bowtie2-build',
             self.reference_fasta_filename(),
             self/'bowtie',
             ],
         stdout = log_to,
         )

def index_vcf(filename):
    """ IGV index a VCF file. 
        Don't fail if igvtools fails (eg not installed).
    """
    try:
        io.execute('igvtools index FILENAME',FILENAME=filename)
    except (OSError, AssertionError):
        print >> sys.stderr, 'Failed to index VCF file with igvtools. Continuing.'

def sort_bam(in_filename, out_prefix, by_name=False, cores=8):
    cores = min(cores, legion.coordinator().get_cores())
    megs = max(10, 800 // cores)
    
    io.execute(
        [ 'samtools', 'sort', '[email protected]', '%d' % cores, '-m', '%dM' % megs ] +
        ([ '-n' ] if by_name else [ ]) +
        [ in_filename, out_prefix ], cores=cores)

    def run(self):
        base = os.path.split(self.prefix)[1]
        
        annotations = [ ]
        sequences = [ ]
        
        for filename in self.filenames:
            any = False
            if io.is_sequence_file(filename):
                sequences.append(filename)
                any = True
            if annotation.is_annotation_file(filename):
                annotations.append(filename)
                any = True
            assert any, 'File is neither a recognized sequence or annotation file'

        cytoband_filename = os.path.join(self.prefix,base+'_cytoband.txt')
        property_filename = os.path.join(self.prefix,'property.txt')
        gff_filename = os.path.join(self.prefix,base+'.gff')
        output_filenames = [ cytoband_filename, property_filename, gff_filename ] 

        if not os.path.exists(self.prefix):
            os.mkdir(self.prefix)
            
        f = open(property_filename,'wb')
        print >> f, 'ordered=true'
        print >> f, 'id=%s' % base
        print >> f, 'name=%s' % (self.name or base)
        print >> f, 'cytobandFile=%s_cytoband.txt' % base
        print >> f, 'geneFile=%s.gff' % base
        print >> f, 'sequenceLocation=%s' % base
        f.close()
        
        trivia.As_gff(output=gff_filename,
               filenames=annotations,
               exclude=[ 'gene', 'source' ]
        ).run()
        
        f_cyt = open(cytoband_filename,'wb')
        for filename in sequences:
            for name, seq in io.read_sequences(filename):
                assert '/' not in name
                f = open(os.path.join(self.prefix, name + '.txt'), 'wb')
                f.write(seq)
                f.close()
                print >> f_cyt, '%s\t0\t%d' % (name, len(seq))
        f_cyt.close()
        
        genome_filename = self.prefix + '.genome'
        if os.path.exists(genome_filename):
            os.unlink(genome_filename)
        io.execute(
            ['zip', '-j', io.abspath(genome_filename)] +
            [ io.abspath(item) for item in output_filenames ]
        )
        for filename in output_filenames:
            if os.path.exists(filename):
                os.unlink(filename)

 def run(self):
     with workspace.tempspace() as temp:
         with open(temp/'batch.txt','wb') as f:
             print >> f, 'new'
             print >> f, 'genome '+os.path.abspath(self.genome)
             for filename in self.files:
                 print >> f, 'load '+os.path.abspath(filename)
         
         io.execute(['java','-jar',io.find_jar('igv.jar'),'-b',temp/'batch.txt'])

 def build_shrimp_mmap(self, cs=False):
     suffix = '-cs' if cs else '-ls'
     
     grace.status('Building SHRiMP mmap')
     io.execute([
         'gmapper' + suffix,
         '--save', self.object_filename('reference' + suffix),
         self.reference_fasta_filename(),
     ])
     grace.status('')

 def run(self):
     from nesoni import io
     
     f_in = self.begin_input()
     f_out = self.begin_output()
     try:
         io.execute(self.command, stdin=f_in, stdout=f_out)
     finally:
         self.end_output(f_out)
         self.end_input(f_in)

 def run(self):
     with workspace.tempspace() as temp:
         with open(temp/'batch.txt','wb') as f:
             print >> f, 'new'
             print >> f, 'preference LAST_TRACK_DIRECTORY', os.getcwd()
             print >> f, 'preference LAST_GENOME_IMPORT_DIRECTORY', os.getcwd()
             print >> f, 'genome '+os.path.abspath(self.genome)
             for filename in self.files:
                 print >> f, 'load '+os.path.abspath(filename)
         
         io.execute(['java','-jar',io.find_jar('igv.jar'),'-b',temp/'batch.txt'])

    def run(self):
        reference = reference_directory.Reference(self.reference, must_exist=True)
        
        jar = io.find_jar('snpEff.jar')
        
        with open(self.prefix + '.vcf','wb') as f:
            io.execute('java -jar JAR eff GENOME VCF -c CONFIG',
                JAR=jar, GENOME=reference.name, VCF=self.vcf, CONFIG=reference/'snpeff.config',
                stdout=f)

        index_vcf(self.prefix+'.vcf')

 def run(self):
     work = self.get_workspace()
     acc = self.run_accession
     io.execute(
         'wget -c URL',
         URL='http://ftp-private.ncbi.nlm.nih.gov/sra/sra-instant/reads/ByRun/sra/%s/%s/%s/%s.sra'
             % (acc[:3],acc[:6],acc,acc),
         cwd=work.working_dir,
         )
     
     io.execute(
         'fastq-dump --split-files --bzip2 FILENAME',
         FILENAME=acc+'.sra',
         cwd=work.working_dir,
         )

 def run(self):
     from nesoni import io
     
     assert self.command, 'Nothing to execute!'
     
     print self.ident()
     
     f_in = self.begin_input()
     f_out = self.begin_output()
     try:
         io.execute(self.command[:1] + self.execution_options + self.command[1:], 
                    stdin=f_in, stdout=f_out)
     finally:
         self.end_output(f_out)
         self.end_input(f_in)

 def run(self):
     assert self.sort in ('queryname', 'coordinate')
     
     jar = io.find_jar('MergeSamFiles.jar', 'MergeSamFiles is part of the Picard package.')
     io.execute([
         'java','-jar',jar,
         'USE_THREADING=true',
         'TMP_DIR='+tempfile.gettempdir(), #Force picard to use same temp dir as Python
         'SORT_ORDER='+self.sort,
         'OUTPUT='+self.prefix+'.bam'
         ] + [ 'INPUT='+item for item in self.bams ])
     
     if self.sort == 'coordinate' and self.index:
         jar = io.find_jar('BuildBamIndex.jar', 'BuildBamIndex is part of the Picard package.')
         io.execute([
             'java','-jar',jar,
             'INPUT='+self.prefix+'.bam'
             ])

 def href(self, filename, title=None, image=False):
     relative = self.workspace.path_as_relative_path(filename)
     if title is None:
        title = os.path.split(filename)[1]
        
        size = os.stat(filename).st_size
        if size >= 1<<30:
            title += ' (%.1fGb)' % (float(size)/(1<<30))
        elif size >= 1<<20: 
            title += ' (%.1fMb)' % (float(size)/(1<<20))
        elif size >= 1<<10: 
            title += ' (%.1fkb)' % (float(size)/(1<<10))
     
     if image:
        thumb_name = 'thumb-'+os.path.splitext(relative)[0]+'.png'
        thumb_filename = self.workspace/thumb_name
        io.execute(['convert', '-thumbnail', '50x50', filename, thumb_filename])
        title = ('<span style="display: inline-block; width: 50px;"><img src="%s"/></span> ' % thumb_name) + title 
        
     return '<a href="%s">%s</a>' % (relative, title)

    def set_sequences(self, filenames):        
        reference_genbank_filename = self / 'reference.gbk'
        reference_filename = self / 'reference.fa'

        reference_genbank_file = open(reference_genbank_filename,'wb')
        any_genbank = [ False ]

        def genbank_callback(name, record):
            """ Make a copy of any genbank files passed in. """
            from Bio import SeqIO
            
            SeqIO.write([record], reference_genbank_file, 'genbank')
            
            f = open(self / (grace.filesystem_friendly_name(name) + '.gbk'), 'wb')
            SeqIO.write([record], f, 'genbank')
            f.close()
            
            any_genbank[0] = True
        
        lengths = [ ]
        seen = set()
        f = open(reference_filename, 'wb')
        for filename in filenames:
            for name, seq in io.read_sequences(filename, genbank_callback=genbank_callback):
                name = name.split()[0]
                assert name not in seen, 'Duplicate chromosome name: ' + name
                seen.add(name)
                lengths.append( (name, len(seq)) )
                io.write_fasta(f, name, seq)
        f.close()        
        self.set_object(lengths, 'reference-lengths.pickle.gz')
        
        reference_genbank_file.close()
        if not any_genbank[0]:
            os.unlink(reference_genbank_filename)
            
        # Create an index of the reference sequences for samtools
        io.execute([
            'samtools', 'faidx', reference_filename
        ])

    def run(self):
        workspace = working_directory.Working(self.output_dir)        
        workspace.setup_reference(self.reference)
        workspace.update_param(snp_cost = self.snp_cost)
        
        #assert os.path.exists(self.reference), 'Reference file does not exist'
        #reference_filename = workspace._object_filename('reference.fa')
        #if os.path.exists(reference_filename):
        #   os.unlink(reference_filename)
        #os.symlink(os.path.relpath(self.reference, self.output_dir), reference_filename)
        
        bam_filename = io.abspath(self.output_dir, 'alignments.bam')
        bam_prefix = io.abspath(self.output_dir, 'alignments')

        if sam.is_bam(self.input):
            sort_input_filename = self.input
            temp_filename = None
        else:
            temp_filename = io.abspath(self.output_dir, 'temp.bam')
            sort_input_filename = temp_filename
            writer = io.Pipe_writer(temp_filename, ['samtools', 'view', '-S', '-b', '-'])
            f = open(self.input, 'rb')
            while True:
                data = f.read(1<<20)
                if not data: break
                writer.write(data)
            writer.close()
            f.close()
        
        grace.status('Sort')
        
        io.execute([
            'samtools', 'sort', '-n', sort_input_filename, bam_prefix
        ])
        
        if temp_filename is not None:
            os.unlink(temp_filename)
        
        grace.status('')

    def run(self):
        workspace = working_directory.Working(self.output_dir)
        workspace.setup_reference(self.reference)

        # assert os.path.exists(self.reference), 'Reference file does not exist'
        # reference_filename = workspace._object_filename('reference.fa')
        # if os.path.exists(reference_filename):
        #   os.unlink(reference_filename)
        # os.symlink(os.path.relpath(self.reference, self.output_dir), reference_filename)

        bam_filename = io.abspath(self.output_dir, "alignments.bam")
        bam_prefix = io.abspath(self.output_dir, "alignments")

        if sam.is_bam(self.input):
            sort_input_filename = self.input
            temp_filename = None
        else:
            temp_filename = io.abspath(self.output_dir, "temp.bam")
            sort_input_filename = temp_filename
            writer = io.Pipe_writer(temp_filename, ["samtools", "view", "-S", "-b", "-"])
            f = open(self.input, "rb")
            while True:
                data = f.read(1 << 20)
                if not data:
                    break
                writer.write(data)
            writer.close()
            f.close()

        grace.status("Sort")

        io.execute(["samtools", "sort", "-n", sort_input_filename, bam_prefix])

        if temp_filename is not None:
            os.unlink(temp_filename)

        grace.status("")