Collect the gaps between adjacent HSPs in the BLAST file.
779017
def collect_gaps(blast, use_subject=False):
    
    key = lambda x: x.sstart if use_subject else x.qstart
    blast.sort(key=key)

    for a, b in zip(blast, blast[1:]):
        if use_subject:
            if a.sstop < b.sstart:
                yield b.sstart - a.sstop
        else:
            if a.qstop < b.qstart:
                yield b.qstart - a.qstop
755912
public static void getAllSplits(final Set<String> splits, final File path, final String trainingSuffix, final String testSuffix) {
        if (path == null) {
            return;
        }
        File[] files = path.listFiles();
        if (files == null) {
            return;
        }
        for (File file : files) {
            if (file.isDirectory()) {
                getAllSplits(splits, file, trainingSuffix, testSuffix);
            } else if (file.getName().endsWith(trainingSuffix)) {
                splits.add(file.getAbsolutePath().replaceAll(trainingSuffix + "$", ""));
            } else if (file.getName().endsWith(testSuffix)) {
                splits.add(file.getAbsolutePath().replaceAll(testSuffix + "$", ""));
            }
        }
    }
959317
def _extract_sections(data_block):
    ''''''
    structure_block = data_block[:, COLS.TYPE:COLS.COL_COUNT].astype(np.int)

    # SWC ID -> structure_block position
    id_map = {-1: -1}
    for i, row in enumerate(structure_block):
        id_map[row[ID]] = i

    # end points have either no children, more than one, or are the start
    # of a new gap
    sec_end_pts = _section_end_points(structure_block, id_map)

    # a 'gap' is when a section has part of it's segments interleaved
    # with those of another section
    gap_sections = set()

    sections = []

    def new_section():
        '''new_section'''
        sections.append(DataBlockSection())
        return sections[-1]

    curr_section = new_section()

    parent_section = {-1: -1}

    for row in structure_block:
        row_id = id_map[row[ID]]
        parent_id = id_map[row[PID]]
        if not curr_section.ids:
            # first in section point is parent
            curr_section.ids.append(parent_id)
            curr_section.ntype = row[TYPE]

        gap = parent_id != curr_section.ids[-1]

        # If parent is not the previous point, create a section end-point.
        # Else add the point to this section
        if gap:
            sec_end_pts.add(row_id)
        else:
            curr_section.ids.append(row_id)

        if row_id in sec_end_pts:
            parent_section[curr_section.ids[-1]] = len(sections) - 1
            # Parent-child discontinuity section
            if gap:
                curr_section = new_section()
                curr_section.ids.extend((parent_id, row_id))
                curr_section.ntype = row[TYPE]
                gap_sections.add(len(sections) - 2)
            elif row_id != len(data_block) - 1:
                # avoid creating an extra DataBlockSection for last row if it's a leaf
                curr_section = new_section()

    for sec in sections:
        # get the section parent ID from the id of the first point.
        if sec.ids:
            sec.pid = parent_section[sec.ids[0]]

        # join gap sections and "disable" first half
        if sec.pid in gap_sections:
            _merge_sections(sections[sec.pid], sec)

    # TODO find a way to remove empty sections.  Currently they are
    # required to maintain tree integrity.
    return sections
1024172
def gaps(self) -> "IntervalList":
        
        if len(self.intervals) < 2:
            return IntervalList(None)
        gaps = []
        for i in range(len(self.intervals) - 1):
            gap = Interval(
                self.intervals[i].end,
                self.intervals[i + 1].start
            )
            gaps.append(gap)
        return IntervalList(gaps)
1327489
def info(args):
    
    p = OptionParser(info.__doc__)
    p.add_option("--gaps", default=False, action="store_true",
                 help="Count number of gaps [default: %default]")
    p.set_table()
    p.set_outfile()
    opts, args = p.parse_args(args)

    if len(args) == 0:
        sys.exit(not p.print_help())

    fastafiles = args
    data = []
    for f in fastafiles:
        s = SequenceInfo(f, gapstats=opts.gaps)
        data.append(s.data)
    write_csv(s.header, data, sep=opts.sep,
              filename=opts.outfile, align=opts.align)
673720
def blastnprep(self):
        """"""
        # Populate threads for each gene, genome combination
        for sample in self.metadata:
            if sample.general.bestassemblyfile != 'NA':
                #
                # sample[self.analysistype].alleleresults = GenObject()
                sample[self.analysistype].closealleles = dict()
                sample[self.analysistype].mismatches = dict()
                sample[self.analysistype].alignmentlength = dict()
                sample[self.analysistype].subjectlength = dict()
                sample[self.analysistype].queryid = dict()
                sample[self.analysistype].start = dict()
                sample[self.analysistype].end = dict()
                sample[self.analysistype].queryseq = dict()
                if type(sample[self.analysistype].allelenames) == list:
                    for allele in sample[self.analysistype].combinedalleles:
                        # Add each fasta/allele file combination to the threads
                        self.runblast(sample.general.bestassemblyfile, allele, sample)
1338154
def gap_index_map(sequence, gap_chars='-'):
    
    return dict(
        (v, k) for k, v in list(ungap_index_map(sequence, gap_chars).items()))
599568
def blast(self, blastfile=None, outfile=None):
        
        from jcvi.formats.blast import BlastSlow, BlastLineByConversion

        if not outfile:
            outfile = self.filename + ".blast"

        if blastfile is not None:
            blasts = BlastSlow(blastfile).to_dict()
        else:
            blasts = None

        fw = must_open(outfile, "w", checkexists=True)
        nlines = 0
        for a, b, id in self.iter_pairs():
            if (a, b) in blasts:
                bline = blasts[(a, b)]
            elif (b, a) in blasts:
                bline = blasts[(b, a)]
            else:
                line = "\t".join((a, b))
                bline = BlastLineByConversion(line, mode="110000000000")

            print(bline, file=fw)
            nlines += 1
        fw.close()

        logging.debug("A total of {0} BLAST lines written to `{1}`."\
                        .format(nlines, outfile))

        return outfile
490619
def get_splits(split_bed, gff_file, stype, key):
    
    bed_file = get_bed_file(gff_file, stype, key)
    cmd = "intersectBed -a {0} -b {1} -wao".format(split_bed, bed_file)
    cmd += " | cut -f4,10"
    p = popen(cmd)
    splits = defaultdict(set)
    for row in p:
        a, b = row.split()
        splits[a].add(b)

    return splits
505867
def read_blast_xml(filename, **kwargs):
    """"""
    # Read file.
    with open(filename, 'r') as f:
        blast_record = NCBIXML.read(f)

    # Prepare DataFrame fields.
    data = {'accession': [],
            'hit_def': [],
            'hit_id': [],
            'title': [],
            'length': [],
            'e_value': [],
            'sequence': []}

    # Get alignments from blast result.
    for i, s in enumerate(blast_record.alignments):
        data['accession'] = s.accession
        data['hit_def'] = s.hit_def
        data['hit_id'] = s.hit_id
        data['title'] = s.title
        data['length'] = s.length
        data['e_value'] = s.hsps[0].expect
        data['sequence'] = s.hsps[0].sbjct

    # Port to DataFrame.
    return pd.DataFrame(data)
215660
def _combine_ranges_on_length(self, data_len, first, second):
        '''
        
        '''
        first = get_true_slice(first, data_len)
        second = get_true_slice(second, data_len)
        final_start, final_step, final_stop = (None, None, None)

        # Get our start
        if first.start == None and second.start == None:
            final_start = None
        else:
            final_start = (first.start if first.start else 0)+(second.start if second.start else 0)

        # Get our stop
        if second.stop == None:
            final_stop = first.stop
        elif first.stop == None:
            final_stop = (first.start if first.start else 0) + second.stop
        else:
            final_stop = min(first.stop, (first.start if first.start else 0) + second.stop)

        # Get our step
        if first.step == None and second.step == None:
            final_step = None
        else:
            final_step = (first.step if first.step else 1)*(second.step if second.step else 1)

        # If we have a start above our stop, set them to be equal
        if final_start > final_stop:
            final_start = final_stop

        return slice(final_start, final_stop, final_step)
1373307
def flanks(args):
    
    p = OptionParser(flanks.__doc__)
    p.add_option("--extend", default=2000, type="int",
                 help="Extend seq flanking the gaps [default: %default]")
    opts, args = p.parse_args(args)

    if len(args) != 2:
        sys.exit(not p.print_help())

    gapsbed, fastafile = args
    Ext = opts.extend
    sizes = Sizes(fastafile).mapping

    bed = Bed(gapsbed)
    pf = gapsbed.rsplit(".", 1)[0]
    extbed = pf + ".ext.bed"
    fw = open(extbed, "w")
    for i, b in enumerate(bed):
        seqid = b.seqid
        gapname = b.accn
        size = sizes[seqid]

        prev_b = bed[i - 1] if i > 0 else None
        next_b = bed[i + 1] if i + 1 < len(bed) else None
        if prev_b and prev_b.seqid != seqid:
            prev_b = None
        if next_b and next_b.seqid != seqid:
            next_b = None

        start = prev_b.end + 1 if prev_b else 1
        start, end = max(start, b.start - Ext), b.start - 1
        print("\t".join(str(x) for x in \
                             (b.seqid, start - 1, end, gapname + "L")), file=fw)

        end = next_b.start - 1 if next_b else size
        start, end = b.end + 1, min(end, b.end + Ext)
        print("\t".join(str(x) for x in \
                             (b.seqid, start - 1, end, gapname + "R")), file=fw)
    fw.close()

    extfasta = fastaFromBed(extbed, fastafile, name=True)
    return extbed, extfasta
561356
def _sufficient_gaps(self,
                         startdate: datetime.date,
                         enddate: datetime.date,
                         requiredgaps: List[datetime.timedelta],
                         flexibility: int) -> Tuple[bool, Optional[Interval]]:
        
        requiredgaps = list(requiredgaps)  # make a copy
        interval = Interval.dayspan(startdate, enddate, include_end=True)
        # log.debug(">>> _sufficient_gaps")
        gaps = self.gap_subset(interval, flexibility)
        gapdurations = gaps.durations()
        gaplist = gaps.list()
        gapdurations.sort(reverse=True)  # longest gap first
        requiredgaps.sort(reverse=True)  # longest gap first
        # log.debug("... gaps = {}".format(gaps))
        # log.debug("... gapdurations = {}".format(gapdurations))
        # log.debug("... requiredgaps = {}".format(requiredgaps))
        while requiredgaps:
            # log.debug("... processing gap")
            if not gapdurations:
                # log.debug("<<< no gaps left")
                return False, None
            if gapdurations[0] < requiredgaps[0]:
                # log.debug("<<< longest gap is too short")
                return False, self.first_interval_ending(gaplist[0].start)
            gapdurations.pop(0)
            requiredgaps.pop(0)
            gaplist.pop(0)
            # ... keeps gaplist and gapdurations mapped to each other
        # log.debug("<<< success")
        return True, None
1346836
private Map<Integer, Integer> getLostStripes(
              Configuration conf, FileStatus stat, FileSystem fs) 
                  throws IOException {
    Map<Integer, Integer> lostStripes = new HashMap<Integer, Integer>();
    RaidInfo raidInfo = RaidUtils.getFileRaidInfo(stat, conf);
    if (raidInfo.codec == null) {
      // Can not find the parity file, the file is not raided.
      return lostStripes;
    }
    Codec codec = raidInfo.codec;
    
    if (codec.isDirRaid) {
      RaidUtils.collectDirectoryCorruptBlocksInStripe(conf, 
          (DistributedFileSystem)fs, raidInfo, 
          stat, lostStripes);
    } else {
      RaidUtils.collectFileCorruptBlocksInStripe((DistributedFileSystem)fs, 
          raidInfo, stat, lostStripes);
    }
    return lostStripes;
  }
746214
private int count_gaps(int[] i1, int[] i2){

		int i0 = i1[0];
		int j0 = i2[0];
		int gaps = 0;
		for (int i =1 ; i<i1.length;i++ ){
			if ( Math.abs(i1[i]-i0) != 1 ||
					( Math.abs(i2[i]-j0) != 1)){
				gaps +=1;
			}
			i0 = i1[i];
			j0 = i2[i];
		}

		return gaps;
	}
1116625
def gap_subset(self, interval: Interval,
                   flexibility: int = 2) -> "IntervalList":
        
        return self.gaps().subset(interval, flexibility)
187044
def separate_globs(globs):
	""""""
	exclude = []
	include = []

	for path in globs:
		if path.startswith("!"):
			exclude.append(path[1:])
		else:
			include.append(path)

	return (exclude, include)