custom = {(reads mapped to scaffold)/(total reads for sample)}/(length of scaffold)
170481
def calc_custom(custom, genome, scaffold, sequence, scaffold_coverage, total_bases):
	
	index = 0
	if scaffold in scaffold_coverage: # what if the scaffold does not have bases mapped back to it? (this *should* not happen)
		if genome not in custom:
			custom[genome] = [[] for i in scaffold_coverage[scaffold]]
		for cov in scaffold_coverage[scaffold]:
			length = float(len(sequence[1]))
			bases = cov * length
			custom_value = ((bases) / (total_bases[index])) / length
			custom[genome][index].append(custom_value)
			index += 1
	return custom
1061582
def read_len( f ):
    """"""
    mapping = dict()
    for line in f:
        fields = line.split()
        mapping[ fields[0] ] = int( fields[1] )
    return mapping
449019
def make_factory(self, cls, count):
        
        field_names = cls._meta.get_all_field_names()
        fields = {}
        text = []
        finalizer = None
        scaffold = scaffolding.scaffold_for_model(cls)

        for field_name in field_names:
            generator = getattr(scaffold, field_name, None)
            if generator:
                if hasattr(generator, 'set_up'):
                    generator.set_up(cls, count)
                fields[field_name] = generator
                text.append(u'%s: %s; ' % (field_name, fields[field_name]))
        try:
            self.stdout.write(u'Generator for %s: %s\n' % (cls, u''.join(text)))
        except models.ObjectDoesNotExist:
            self.stdout.write(u'Generator for %s\n' % u''.join(text))

        if hasattr(scaffold, 'finalize') and hasattr(scaffold.finalize, '__call__'):
            finalizer = scaffold.finalize

        return fields, finalizer
749224
def split_input(cls, mapper_spec):
    
    params = _get_params(mapper_spec)
    shard_count = mapper_spec.shard_count

    # Pick out the overall start and end times and time step per shard.
    start_time = params[cls.START_TIME_PARAM]
    end_time = params[cls.END_TIME_PARAM]
    seconds_per_shard = (end_time - start_time) / shard_count

    # Create a LogInputReader for each shard, modulating the params as we go.
    shards = []
    for _ in xrange(shard_count - 1):
      params[cls.END_TIME_PARAM] = (params[cls.START_TIME_PARAM] +
                                    seconds_per_shard)
      shards.append(LogInputReader(**params))
      params[cls.START_TIME_PARAM] = params[cls.END_TIME_PARAM]

    # Create a final shard to complete the time range.
    params[cls.END_TIME_PARAM] = end_time
    return shards + [LogInputReader(**params)]
324647
def main():
  
  r = Random(42)

  # Measure serialization
  startSerializationTime = time.time()

  for i in xrange(_SERIALIZATION_LOOPS):
    # NOTE pycapnp's builder.from_dict (used in nupic.bindings) leaks
    # memory if called on the same builder more than once, so we construct a
    # fresh builder here
    builderProto = RandomProto.new_message()
    r.write(builderProto)

  elapsedSerializationTime = time.time() - startSerializationTime


  builderBytes = builderProto.to_bytes()


  # Measure deserialization
  startDeserializationTime = time.time()

  deserializationCount = 0
  while deserializationCount < _DESERIALIZATION_LOOPS:
    # NOTE: periodicaly create a new reader to avoid "Exceeded message traversal
    # limit" error
    readerProto = RandomProto.from_bytes(
      builderBytes,
      traversal_limit_in_words=_TRAVERSAL_LIMIT_IN_WORDS,
      nesting_limit=_NESTING_LIMIT)

    numReads = min(_DESERIALIZATION_LOOPS - deserializationCount,
                     _MAX_DESERIALIZATION_LOOPS_PER_READER)
    for _ in xrange(numReads):
      r.read(readerProto)

    deserializationCount += numReads

  elapsedDeserializationTime = time.time() - startDeserializationTime

  # Print report
  print _SERIALIZATION_LOOPS, "Serialization loops in", \
        elapsedSerializationTime, "seconds."
  print "\t", elapsedSerializationTime/_SERIALIZATION_LOOPS, "seconds per loop."

  print deserializationCount, "Deserialization loops in", \
        elapsedDeserializationTime, "seconds."
  print "\t", elapsedDeserializationTime/deserializationCount, "seconds per loop."
617540
private LocalVariableDefinition constructStructInstance(MethodDefinition read, Map<Short, LocalVariableDefinition> structData)
    {
        LocalVariableDefinition instance = read.addLocalVariable(structType, "instance");

        // create the new instance (or builder)
        if (metadata.getBuilderClass() == null) {
            read.newObject(structType).dup();
        }
        else {
            read.newObject(metadata.getBuilderClass()).dup();
        }

        // invoke constructor
        ThriftConstructorInjection constructor = metadata.getConstructorInjection().get();
        // push parameters on stack
        for (ThriftParameterInjection parameter : constructor.getParameters()) {
            read.loadVariable(structData.get(parameter.getId()));
        }
        // invoke constructor
        read.invokeConstructor(constructor.getConstructor())
                .storeVariable(instance);
        return instance;
    }
721012
def simulate_reads(self):
        
        logging.info('Read simulation')
        for sample in self.metadata:
            # Create the simulated_reads GenObject
            sample.simulated_reads = GenObject()
            # Iterate through all the desired depths of coverage
            for depth in self.read_depths:
                # Create the depth GenObject
                setattr(sample.simulated_reads, depth, GenObject())
                # Set the depth and output directory attributes for the depth GenObject
                sample.simulated_reads[depth].depth = depth
                sample.simulated_reads[depth].depth_dir = os.path.join(sample.outputdir, 'simulated', depth)
                # Create the output directory
                make_path(sample.simulated_reads[depth].depth_dir)
                # Iterate through all the desired forward and reverse read pair lengths
                for read_pair in self.read_lengths:
                    # Create the read_pair GenObject within the depth GenObject
                    setattr(sample.simulated_reads[depth], read_pair, GenObject())
                    # Set and create the output directory
                    sample.simulated_reads[depth][read_pair].outputdir = \
                        os.path.join(sample.simulated_reads[depth].depth_dir, read_pair)
                    make_path(sample.simulated_reads[depth][read_pair].outputdir)
                    # Create both forward_reads and reverse_reads sub-GenObjects
                    sample.simulated_reads[depth][read_pair].forward_reads = GenObject()
                    sample.simulated_reads[depth][read_pair].reverse_reads = GenObject()
                    # Extract the forward and reverse reads lengths from the read_pair variable
                    sample.simulated_reads[depth][read_pair].forward_reads.length, \
                        sample.simulated_reads[depth][read_pair].reverse_reads.length = read_pair.split('_')
                    # Set the name of the forward reads - include the depth and read length information
                    sample.simulated_reads[depth][read_pair].forward_reads.fastq = \
                        os.path.join(sample.simulated_reads[depth][read_pair].outputdir,
                                     '{name}_{depth}_{read_pair}_R1.fastq.gz'
                                     .format(name=sample.name,
                                             depth=depth,
                                             read_pair=read_pair))
                    # Reverse reads
                    sample.simulated_reads[depth][read_pair].reverse_reads.fastq = \
                        os.path.join(sample.simulated_reads[depth][read_pair].outputdir,
                                     '{name}_{depth}_{read_pair}_R2.fastq.gz'
                                     .format(name=sample.name,
                                             depth=depth,
                                             read_pair=read_pair))
                    # Create the trimmed output directory attribute
                    sample.simulated_reads[depth][read_pair].simulated_trimmed_outputdir \
                        = os.path.join(sample.simulated_reads[depth][read_pair].outputdir,
                                       'simulated_trimmed')
                    # Set the name of the forward trimmed reads - include the depth and read length information
                    # This is set now, as the untrimmed files will be removed, and a check is necessary
                    sample.simulated_reads[depth][read_pair].forward_reads.trimmed_simulated_fastq = \
                        os.path.join(sample.simulated_reads[depth][read_pair].simulated_trimmed_outputdir,
                                     '{name}_simulated_{depth}_{read_pair}_R1.fastq.gz'
                                     .format(name=sample.name,
                                             depth=depth,
                                             read_pair=read_pair))
                    # Reverse reads
                    sample.simulated_reads[depth][read_pair].reverse_reads.trimmed_simulated_fastq = \
                        os.path.join(sample.simulated_reads[depth][read_pair].simulated_trimmed_outputdir,
                                     '{name}_simulated_{depth}_{read_pair}_R2.fastq.gz'
                                     .format(name=sample.name,
                                             depth=depth,
                                             read_pair=read_pair))
                    # Calculate the number of reads required for the forward and reverse reads to yield the
                    # desired coverage depth e.g. 5Mbp genome at 20X coverage: 100Mbp in reads. 50bp forward reads
                    # 150bp reverse reads: forward proportion is 50 / (150 + 50) = 0.25 (and reverse is 0.75).
                    # Forward total reads is 25Mbp (75Mbp reverse). Number of reads required = 25Mbp / 50 bp
                    # 500000 reads total (same for reverse, as the reads are longer)
                    sample.simulated_reads[depth][read_pair].num_reads = \
                        int(sample.assembly_length *
                            int(depth) *
                            (int(sample.simulated_reads[depth][read_pair].forward_reads.length) /
                             (int(sample.simulated_reads[depth][read_pair].forward_reads.length) +
                              int(sample.simulated_reads[depth][read_pair].reverse_reads.length)
                              )
                             ) /
                            int(sample.simulated_reads[depth][read_pair].forward_reads.length)
                            )
                    logging.info(
                        'Simulating {num_reads} paired reads for sample {name} with the following parameters:\n'
                        'depth {dp}, forward reads {fl}bp, and reverse reads {rl}bp'
                        .format(num_reads=sample.simulated_reads[depth][read_pair].num_reads,
                                dp=depth,
                                name=sample.name,
                                fl=sample.simulated_reads[depth][read_pair].forward_reads.length,
                                rl=sample.simulated_reads[depth][read_pair].reverse_reads.length))
                    # If the reverse reads are set to 0, supply different parameters to randomreads
                    if sample.simulated_reads[depth][read_pair].reverse_reads.length != '0':
                        # Ensure that both the simulated reads, and the trimmed simulated reads files don't
                        # exist before simulating the reads
                        if not os.path.isfile(sample.simulated_reads[depth][read_pair].forward_reads.fastq) and \
                                not os.path.isfile(
                                    sample.simulated_reads[depth][read_pair].forward_reads.trimmed_simulated_fastq):
                            # Use the randomreads method in the OLCTools bbtools wrapper to simulate the reads
                            out, \
                                err, \
                                sample.simulated_reads[depth][read_pair].forward_reads.simulate_call = bbtools\
                                .randomreads(reference=sample.bestassemblyfile,
                                             length=sample.simulated_reads[depth][read_pair].reverse_reads.length,
                                             reads=sample.simulated_reads[depth][read_pair].num_reads,
                                             out_fastq=sample.simulated_reads[depth][read_pair].forward_reads.fastq,
                                             paired=True,
                                             returncmd=True,
                                             **{'ziplevel': '9',
                                                'illuminanames': 't',
                                                'Xmx': self.mem}
                                             )
                        else:
                            try:
                                forward_size = os.path.getsize(sample.simulated_reads[depth][read_pair]
                                                               .forward_reads.fastq)
                            except FileNotFoundError:
                                forward_size = 0
                            try:
                                reverse_size = os.path.getsize(sample.simulated_reads[depth][read_pair]
                                                               .reverse_reads.fastq)
                            except FileNotFoundError:
                                reverse_size = 0
                            if forward_size <= 100 or reverse_size <= 100:
                                try:
                                    os.remove(sample.simulated_reads[depth][read_pair].forward_reads.fastq)
                                except FileNotFoundError:
                                    pass
                                try:
                                    os.remove(sample.simulated_reads[depth][read_pair].reverse_reads.fastq)
                                except FileNotFoundError:
                                    pass
                                # Use the randomreads method in the OLCTools bbtools wrapper to simulate the reads
                                out, \
                                err, \
                                sample.simulated_reads[depth][read_pair].forward_reads.simulate_call = bbtools \
                                    .randomreads(reference=sample.bestassemblyfile,
                                                 length=sample.simulated_reads[depth][read_pair].reverse_reads.length,
                                                 reads=sample.simulated_reads[depth][read_pair].num_reads,
                                                 out_fastq=sample.simulated_reads[depth][read_pair].forward_reads.fastq,
                                                 paired=True,
                                                 returncmd=True,
                                                 **{'ziplevel': '9',
                                                    'illuminanames': 't'}
                                                 )
                    else:
                        if not os.path.isfile(sample.simulated_reads[depth][read_pair].forward_reads.fastq):
                            # Use the randomreads method in the OLCTools bbtools wrapper to simulate the reads
                            out, \
                                err, \
                                sample.simulated_reads[depth][read_pair].forward_reads.simulate_call = bbtools\
                                .randomreads(reference=sample.bestassemblyfile,
                                             length=sample.simulated_reads[depth][read_pair].forward_reads.length,
                                             reads=sample.simulated_reads[depth][read_pair].num_reads,
                                             out_fastq=sample.simulated_reads[depth][read_pair].forward_reads.fastq,
                                             paired=False,
                                             returncmd=True,
                                             **{'ziplevel': '9',
                                                'illuminanames': 't'}
                                             )
                # Update the JSON file
                self.write_json(sample)
1177125
def _LinearFoldByteStream(self, mapped_value, **unused_kwargs):
    
    try:
      attribute_values = [
          getattr(mapped_value, attribute_name, None)
          for attribute_name in self._attribute_names]
      attribute_values = [
          value for value in attribute_values if value is not None]
      return self._operation.WriteTo(tuple(attribute_values))

    except Exception as exception:
      error_string = (
          'Unable to write: {0:s} to byte stream with error: {1!s}').format(
              self._data_type_definition.name, exception)
      raise errors.FoldingError(error_string)
370270
def Read(self, length=None):
    """"""

    if self.progress_callback:
      self.progress_callback()

    available_to_read = max(0, (self.size or 0) - self.offset)
    if length is None:
      to_read = available_to_read
    else:
      to_read = min(length, available_to_read)

    with FileHandleManager(self.filename) as fd:
      offset = self.file_offset + self.offset
      pre_padding = offset % self.alignment

      # Due to alignment we read some more data than we need to.
      aligned_offset = offset - pre_padding

      fd.Seek(aligned_offset)

      data = fd.Read(to_read + pre_padding)
      self.offset += len(data) - pre_padding

      return data[pre_padding:]
1232888
def sample(self, n):
    

    row_total_count = 0
    row_counts = []
    for file in self.files:
      with _util.open_local_or_gcs(file, 'r') as f:
        num_lines = sum(1 for line in f)
        row_total_count += num_lines
        row_counts.append(num_lines)

    names = None
    dtype = None
    if self._schema:
      _MAPPINGS = {
        'FLOAT': np.float64,
        'INTEGER': np.int64,
        'TIMESTAMP': np.datetime64,
        'BOOLEAN': np.bool,
      }
      names = [x['name'] for x in self._schema]
      dtype = {x['name']: _MAPPINGS.get(x['type'], object) for x in self._schema}

    skip_count = row_total_count - n
    # Get all skipped indexes. These will be distributed into each file.
    # Note that random.sample will raise Exception if skip_count is greater than rows count.
    skip_all = sorted(random.sample(range(0, row_total_count), skip_count))
    dfs = []
    for file, row_count in zip(self.files, row_counts):
      skip = [x for x in skip_all if x < row_count]
      skip_all = [x - row_count for x in skip_all if x >= row_count]
      with _util.open_local_or_gcs(file, 'r') as f:
        dfs.append(pd.read_csv(f, skiprows=skip, names=names, dtype=dtype, header=None))
    return pd.concat(dfs, axis=0, ignore_index=True)
1044224
def parse_cov(cov_table, scaffold2genome):
    
    size   = {} # size[genome] = genome size
    mapped = {} # mapped[genome][sample] = mapped bases
    # parse coverage files
    for line in open(cov_table):
        line = line.strip().split('\t')
        if line[0].startswith('#'):
            samples = line[1:]
            samples = [i.rsplit('/', 1)[-1].split('.', 1)[0] for i in samples]
            continue
        scaffold, length = line[0].split(': ')
        length = float(length)
        covs  = [float(i) for i in line[1:]]
        bases = [c * length for c in covs]
        if scaffold not in scaffold2genome:
            continue
        genome = scaffold2genome[scaffold]
        if genome not in size:
            size[genome] = 0
            mapped[genome] = {sample:0 for sample in samples}
        # keep track of genome size
        size[genome] += length
        # keep track of number of mapped bases
        for sample, count in zip(samples, bases):
            mapped[genome][sample] += count
    # calculate coverage from base counts and genome size
    coverage = {'genome':[], 'genome size (bp)':[], 'sample':[], 'coverage':[]}
    for genome, length in size.items():
        for sample in samples:
            cov = mapped[genome][sample] / length
            coverage['genome'].append(genome)
            coverage['genome size (bp)'].append(length)
            coverage['sample'].append(sample)
            coverage['coverage'].append(cov)
    return pd.DataFrame(coverage)
775719
def _ReadEncodedData(self, read_size):
    
    encoded_data = self._file_object.read(read_size)

    read_count = len(encoded_data)

    self._encoded_data = b''.join([self._encoded_data, encoded_data])

    self._decoded_data, self._encoded_data = (
        self._decoder.Decode(self._encoded_data))

    self._decoded_data_size = len(self._decoded_data)

    return read_count
873833
def _reads_per_position(bam_in, loci_file, out_dir):
    
    data = Counter()
    a = pybedtools.BedTool(bam_in)
    b = pybedtools.BedTool(loci_file)
    c = a.intersect(b, s=True, bed=True, wo=True)
    for line in c:
        end = int(line[1]) + 1 + int(line[2]) if line[5] == "+" else int(line[1]) + 1
        start = int(line[1]) + 1 if line[5] == "+" else int(line[1]) + 1 + int(line[2])
        side5 = "%s\t5p\t%s" % (line[15], start)
        side3 = "%s\t3p\t%s" % (line[15], end)
        data[side5] += 1
        data[side3] += 1

    counts_reads = op.join(out_dir, 'locus_readpos.counts')
    with open(counts_reads, 'w') as out_handle:
        for k in data:
            print(k, file=out_handle, end="")

    return counts_reads
4955
def copies(mapping, s2bins, rna, min_rna = 800, mismatches = 0):
    
    cov = {} # cov[scaffold] = [bases, length]
    s2bins, bins2s = parse_s2bins(s2bins)
    rna_cov = parse_rna(rna, s2bins, min_rna)
    s2bins, bins2s = filter_missing_rna(s2bins, bins2s, rna_cov)
    # count bases mapped to scaffolds and rRNA gene regions
    for line in mapping:
        line = line.strip().split()
        # get scaffold lengths
        if line[0].startswith('@'):
            if line[0].startswith('@SQ') is False:
                continue
            s = line[1].split(':')[1]
            l = int(line[2].split(':')[1])
            # check if scaffold is binned
            if s not in s2bins:
                continue
            if s not in cov:
                cov[s] = [0, l]
        # check mismatch threshold
        mm = count_mismatches(line)
        if mm is False or mm > mismatches:
            continue
        # check that scaffold is in bin
        s, bases = line[2], len(line[9])
        if s not in cov:
            continue
        cov[s][0] += bases 
        rna_cov = rna_bases(rna_cov, s, bases, line) 
    print('# mismatches threshold: %s' % (mismatches))
    header = ['#rRNA scaffold', 'rRNA genes >=%sbp on scaffold' % (min_rna), \
            'rRNA coverage', \
            'bin', 'bin info', 'bin coverage', \
            'rRNAs >=%sbp in bin' % (min_rna), \
            'rRNA coverage/bin coverage', \
            'estimated number of copies']
    print('\t'.join(header))
    for bin, scaffolds in list(bins2s.items()):
        rna_count = sum([len(rna_cov[s][2]) for s in scaffolds if s in rna_cov])
        for s in scaffolds:
            if s not in rna_cov:
                continue
            out = []
            counts = rna_cov[s]
            bin_cov = calc_bin_cov(bins2s[bin], cov)
            num_genes = len(counts[2])
            rna_coverage = float(float(counts[0])/float(counts[1]))
            if bin_cov == 0:
                rna_div_bin = 0
            else:
                rna_div_bin = float(rna_coverage/bin_cov)
            est = int(max([rna_count, counts, rna_div_bin]))
            out = [s, num_genes, rna_coverage, bin, bin_cov, rna_count, rna_div_bin, est]
            print('\t'.join([str(i) for i in out]))
630350
static OrcParseSetup deriveParseSetup(Reader orcFileReader, StructObjectInspector insp) {
    List<StructField> allColumns = (List<StructField>) insp.getAllStructFieldRefs();  // grab column info
    List<StripeInformation> allStripes = orcFileReader.getStripes();  // grab stripe information
    ArrayList<String> allColNames = new ArrayList<>();
    boolean[] toInclude = new boolean[allColumns.size()+1];
    int supportedFieldCnt = 0 ;
    int colIdx = 0;
    for (StructField oneField:allColumns) {
      allColNames.add(oneField.getFieldName());
      String columnType = oneField.getFieldObjectInspector().getTypeName();
      if (columnType.toLowerCase().contains("decimal")) {
        columnType = "decimal";
      }
      if (isSupportedSchema(columnType)) {
        toInclude[colIdx+1] = true;
        supportedFieldCnt++;
      }
      int cnt = countStructFields(oneField.getFieldObjectInspector(),allColNames);
      if(cnt > 1)
        toInclude = Arrays.copyOf(toInclude,toInclude.length + cnt-1);
      colIdx+=cnt;
    }
    String [] allNames = allColNames.toArray(new String[allColNames.size()]);
    String[] names = new String[supportedFieldCnt];

    byte[] types = new byte[supportedFieldCnt];
    String[][] domains = new String[supportedFieldCnt][];
    String[] dataPreview = new String[supportedFieldCnt];
    String[] dataTypes = new String[supportedFieldCnt];
    ParseWriter.ParseErr[] errs = new ParseWriter.ParseErr[0];

    // go through all column information
    int columnIndex = 0;
    for (StructField oneField : allColumns) {
      String columnType = oneField.getFieldObjectInspector().getTypeName();
      if (columnType.toLowerCase().contains("decimal"))
        columnType = "decimal"; // get rid of strange attachment
      if (isSupportedSchema(columnType)) {
        names[columnIndex] = oneField.getFieldName();
        types[columnIndex] = schemaToColumnType(columnType);
        dataTypes[columnIndex] = columnType;
        columnIndex++;
      } else {
        errs = ArrayUtils.append(errs, new ParseWriter.ParseErr("Orc Parser: Skipping field: "
            + oneField.getFieldName() + " because of unsupported type: " + columnType, -1, -1L, -2L));
      }
    }

    // get size of each stripe
    long[] stripeSizes = new long[allStripes.size()];
    long fileSize = 0L;
    long maxStripeSize = 0L;

    for (int index = 0; index < allStripes.size(); index++) {
      long stripeSize = allStripes.get(index).getDataLength();

      if (stripeSize > maxStripeSize)
        maxStripeSize = stripeSize;

      fileSize = fileSize + stripeSize;
      stripeSizes[index] = fileSize;
    }
    OrcParseSetup ps = new OrcParseSetup(
        supportedFieldCnt,
        names,
        types,
        domains,
        null,
        new String[][] { dataPreview },
        orcFileReader,
        dataTypes,
        toInclude,
        allNames,
        errs
    );

    return ps;
  }
883897
def __shpFileLength(self):
        """"""
        # Start with header length
        size = 100
        # Calculate size of all shapes
        for s in self._shapes:
            # Add in record header and shape type fields
            size += 12
            # nParts and nPoints do not apply to all shapes
            #if self.shapeType not in (0,1):
            #       nParts = len(s.parts)
            #       nPoints = len(s.points)
            if hasattr(s,'parts'):
                nParts = len(s.parts)
            if hasattr(s,'points'):
                nPoints = len(s.points)
            # All shape types capable of having a bounding box
            if self.shapeType in (3,5,8,13,15,18,23,25,28,31):
                size += 32
            # Shape types with parts
            if self.shapeType in (3,5,13,15,23,25,31):
                # Parts count
                size += 4
                # Parts index array
                size += nParts * 4
            # Shape types with points
            if self.shapeType in (3,5,8,13,15,23,25,31):
                # Points count
                size += 4
                # Points array
                size += 16 * nPoints
            # Calc size of part types for Multipatch (31)
            if self.shapeType == 31:
                size += nParts * 4
            # Calc z extremes and values
            if self.shapeType in (13,15,18,31):
                # z extremes
                size += 16
                # z array
                size += 8 * nPoints
            # Calc m extremes and values
            if self.shapeType in (23,25,31):
                # m extremes
                size += 16
                # m array
                size += 8 * nPoints
            # Calc a single point
            if self.shapeType in (1,11,21):
                size += 16
            # Calc a single Z value
            if self.shapeType == 11:
                size += 8
            # Calc a single M value
            if self.shapeType in (11,21):
                size += 8
        # Calculate size as 16-bit words
        size //= 2
        return size
384853
def calc_bin_cov(scaffolds, cov):
    
    bases = sum([cov[i][0] for i in scaffolds if i in cov])
    length = sum([cov[i][1] for i in scaffolds if i in cov])
    if length == 0:
        return 0
    return float(float(bases)/float(length))