Java Code Examples for htsjdk.samtools.SAMRecord#getReferenceName()

/**
 * For each alignment block, see which genes have exons that intersect it.
 * Only count genes where the alignment blocks are consistent - either the alignment blocks point to the same gene, or
 * an alignment block points to a gene and the other to no genes/exons in the set.
 * Note that this is not strand specific, which may cause problems if a read overlaps two genes on opposite strands that have overlapping exons.
 * @param rec
 * @param genes A set of genes that the read originally overlaps.
 * @param allowMultipleGenes.  If false, and a read overlaps multiple gene exons, then none of the genes are returned.  If true, return the set of all genes.
 * @return
 */
//TODO: if this is used in the future, make it strand specific.
public Set<Gene> getConsistentExons (final SAMRecord rec, final Set<Gene> genes, final boolean allowMultiGeneReads) {

	Set<Gene> result = new HashSet<>();
	String refName = rec.getReferenceName();
	List<AlignmentBlock> alignmentBlocks = rec.getAlignmentBlocks();
	for (AlignmentBlock b: alignmentBlocks) {
		Set<Gene> blockGenes = getAlignmentBlockonGeneExon(refName, b, genes);
		result.addAll(blockGenes);
		/*
		// if result is not empty and blockGenes isn't empty, intersect the set and set the result as this new set.
		if (result.size()>0 && blockGenes.size()>0) {
			if (allowMultiGeneReads) result.addAll(blockGenes);
			if (!allowMultiGeneReads)  result.retainAll(blockGenes);
		} else
			// if blockGenes is populated and you're here, then result is empty, so set result to these results
			result=blockGenes;
		*/
	}
	if (!allowMultiGeneReads & result.size()>1) return new HashSet<>();
	return result;
}
 

private SAMRecord tagRead(final SAMRecord r, final OverlapDetector<Interval> od) {
	Set<Interval>intervals = new HashSet<>();
	List<AlignmentBlock> blocks = r.getAlignmentBlocks();
	for (AlignmentBlock b: blocks) {
		int refStart =b.getReferenceStart();
		int refEnd = refStart+b.getLength()-1;
		Interval v = new Interval(r.getReferenceName(), refStart, refEnd);
		Collection<Interval> blockResult = od.getOverlaps(v);
		intervals.addAll(blockResult);
	}
	String tagName = getIntervalName(intervals);
	if (tagName!=null)
		r.setAttribute(this.TAG, tagName);
	else
		r.setAttribute(this.TAG, null);
	return (r);

}
 

/**
 * Check if a read overlaps any SNPs in the OverlapDetector.  Tag reads with SNPs.
 * If more than 1 SNP tags a read, make a read for each SNP.
 * Simplified since data goes through GeneFunctionIteratorWrapper to take care of how reads/genes interact.
 */
private void processSNP (final SAMRecord r) {
	List<AlignmentBlock> blocks = r.getAlignmentBlocks();

	Collection<String> snps = new HashSet<>();

	for (AlignmentBlock b: blocks) {
		int start = b.getReferenceStart();
		int end = start + b.getLength() -1;

		Interval i = new Interval(r.getReferenceName(), start, end);
		Collection<Interval> overlaps = this.snpIntervals.getOverlaps(i);
		for (Interval o: overlaps)
			snps.add(IntervalTagComparator.toString(o));
	}

	// 1 read per SNP.
	for (String snp:snps) {
		SAMRecord rr = Utils.getClone(r);
		rr.setAttribute(this.snpTag, snp);
		queueRecordForOutput(rr);
	}
}
 

public Map<Gene, LocusFunction> getLocusFunctionForReadByGene (final SAMRecord rec, final OverlapDetector<Gene> geneOverlapDetector) {
	Map<Gene, LocusFunction> result = new HashMap<>();
	final Interval readInterval = new Interval(rec.getReferenceName(), rec.getAlignmentStart(), rec.getAlignmentEnd(), rec.getReadNegativeStrandFlag(), rec.getReadName());
	final Collection<Gene> overlappingGenes = geneOverlapDetector.getOverlaps(readInterval);

       for (Gene g: overlappingGenes) {
       	LocusFunction f = getLocusFunctionForRead(rec, g);
           result.put(g, f);
       }
       return result;
}
 

public Object getValue(SAMRecord record) {
	switch (type) {
	case QNAME:
		return record.getReadName();
	case FLAG:
		return Integer.toString(record.getFlags());
	case RNAME:
		return record.getReferenceName();
	case POS:
		return Integer.toString(record.getAlignmentStart());
	case MAPQ:
		return Integer.toString(record.getMappingQuality());
	case CIGAR:
		return record.getCigarString();
	case RNEXT:
		return record.getMateReferenceName();
	case PNEXT:
		return Integer.toString(record.getMateAlignmentStart());
	case TLEN:
		return Integer.toString(record.getInferredInsertSize());
	case SEQ:
		return record.getReadString();
	case QUAL:
		return record.getBaseQualityString();

	case TAG:
		if (tagId == null)
			throw new IllegalArgumentException("Tag mismatch reqiues tag id. ");
		return record.getAttribute(tagId);

	default:
		throw new IllegalArgumentException("Unknown record field: " + type.name());
	}
}
 

public static Object getValue(SAMRecord record, FIELD_TYPE field, String tagId) {
	if (field == null)
		throw new IllegalArgumentException("Record field is null.");

	switch (field) {
	case QNAME:
		return record.getReadName();
	case FLAG:
		return Integer.toString(record.getFlags());
	case RNAME:
		return record.getReferenceName();
	case POS:
		return Integer.toString(record.getAlignmentStart());
	case MAPQ:
		return Integer.toString(record.getMappingQuality());
	case CIGAR:
		return record.getCigarString();
	case RNEXT:
		return record.getMateReferenceName();
	case PNEXT:
		return Integer.toString(record.getMateAlignmentStart());
	case TLEN:
		return Integer.toString(record.getInferredInsertSize());
	case SEQ:
		return record.getReadString();
	case QUAL:
		return record.getBaseQualityString();

	case TAG:
		if (tagId == null)
			throw new IllegalArgumentException("Tag mismatch reqiues tag id. ");
		return record.getAttribute(tagId);

	default:
		throw new IllegalArgumentException("Unknown record field: " + field.name());
	}
}
 

SameCoordReadSet getSameCoordReads(PeekIterator<SAMRecord> it, String fileName) throws Exception {
  SameCoordReadSet ret = null;
  try {
    SAMRecord record;
    while (it.hasNext()) {
      record = it.peek();
      if (record.isSecondaryOrSupplementary() ||
          (options.ignoreUnmappedReads && record.getReadUnmappedFlag())) {
        it.next();
        continue;
      }
      if (ret != null) {
        if (record.getAlignmentStart() != ret.coord || !record.getReferenceName().equals(ret.reference)) {
          break;
        }
      } else {
        ret = new SameCoordReadSet();
        ret.map = Maps.newHashMap();
        ret.coord = record.getAlignmentStart();
        ret.reference = record.getReferenceName();
      }
      ret.map.put(record.getReadName(), record);
      it.next();
    }
  } catch (Exception ex) {
    throw new Exception("Error reading from " + fileName + "\n" + ex.getMessage());
  }
  return ret;
}
 

@Override
public SAMRecord next() {

	SAMRecord read = null;
	boolean isCacheUpToDate = false;
	if (!cache.isEmpty()) {
		InsertShiftSAMRecord first = cache.first();
		InsertShiftSAMRecord last = cache.last();
		
		// Don't seek too far ahead
		if (last.read.getAlignmentStart() > first.read.getAlignmentStart()+2 || !last.read.getReferenceName().equals(first.read.getReferenceName())) {
			isCacheUpToDate = true;
		}
		
		read = first.read;
	} else {
		read = iter.next();
		cache.add(new InsertShiftSAMRecord(read));
	}
	
	int cacheStart = read.getAlignmentStart() + 1;
	String cacheChromosome = read.getReferenceName();
	
	while (!isCacheUpToDate && iter.hasNext() && read.getAlignmentStart() <= cacheStart+1 && read.getReferenceName().equals(cacheChromosome)) {
		read = iter.next();
		cache.add(new InsertShiftSAMRecord(read));
	}
	
	return cache.pollFirst().read;
}
 

private boolean getCachedReadsAtCurrentLocus(List<SAMRecord> reads) {
	
	reads.clear();
	
	Iterator<SAMRecord> cacheIter = readCache.iterator();
	
	String nextChr = null;
	int nextPos = -1;
	
	while (cacheIter.hasNext()) {
		SAMRecord read = cacheIter.next();
		
		if (read.getAlignmentEnd() < currentPos && read.getReferenceName().equals(currentChr)) {
			// We've gone past the end of this read, so remove from cache.
			cacheIter.remove();
		} else if (getAlignmentStart(read) <= currentPos && read.getAlignmentEnd() >= currentPos) {
			// This read spans the current locus of interest.
			reads.add(read);
		} else {
			// This read is beyond the current locus.
			if (nextChr == null) {
				nextChr = read.getReferenceName();
				nextPos = getAlignmentStart(read);
			}
		}
	}
	
	if (reads.isEmpty() && nextChr != null) {
		currentChr = nextChr;
		currentPos = nextPos;
		
		return true;
	}
	
	return false;
}
 

private void loadReadsIntoCache() {
	boolean shouldReadFromFile = false;
	
	if (readCache.isEmpty()) {
		shouldReadFromFile = true;
	}
	else {
		SAMRecord last = readCache.get(readCache.size()-1);
		if (getAlignmentStart(last) <= currentPos && last.getReferenceName().equals(currentChr)) {
			shouldReadFromFile = true;
		}
	}
	
	while (shouldReadFromFile && samIter.hasNext()) {
		SAMRecord read = samIter.next();

		// Skip over unmapped reads
		if (!read.getReadUnmappedFlag()) {
		
			if (readCache.isEmpty() && !read.getReferenceName().equals(currentChr)) {
				currentChr = read.getReferenceName();
				currentPos = getAlignmentStart(read);
			}
			
			readCache.add(read);
			
			if (getAlignmentStart(read) > currentPos || !read.getReferenceName().equals(currentChr)) {
				shouldReadFromFile = false;
			}
		}
	}
	
	// Skip huge pileups!
	if (readCache.size() > maxDepth) {
		Logger.warn("Depth too high, clearing read cache " + currentChr + ":" + currentPos);
		for (int i=readCache.size()-2; i>=0; i--) {
			readCache.remove(i);
		}
	}
}
 

public static ReadRecord from(final SAMRecord record)
{
    ReadRecord read = new ReadRecord(
            record.getReadName(), record.getReferenceName(), record.getStart(), record.getEnd(),
            record.getReadString(), record.getCigar(), record.getInferredInsertSize(), record.getFlags(),
            record.getMateReferenceName(), record.getMateAlignmentStart());

    read.setSuppAlignment(record.getStringAttribute(SUPPLEMENTARY_ATTRIBUTE));
    return read;
}
 

/**
 * Does the reference have a soft match for any of the proposed matches?
 * @param matches the list of possible matches
 * @param r the read
 * @return return true if the record matches any of the filters.
 */
boolean exactMatchReference (final List<String> matches, final SAMRecord r) {
	String refName = r.getReferenceName();
	for (String match: matches)
		if (refName.matches(match))
			return true;
	return (false);

}
 

/**
 * Does the reference have a soft match for any of the proposed matches?
 * @param matches the list of possible matches.  All matches are wrapped with .* on either side.
 * @param r the read
 * @return return true if the record matches any of the filters.
 */
boolean softMatchReference (final List<String> matches, final SAMRecord r) {
	String refName = r.getReferenceName();
	for (String match: matches)
		if (refName.matches(".*"+match+".*"))
			return true;
	return (false);

}
 

@Override
    public SAMRecord next() {
        final SAMRecord rec = it.next();
        final String reference = rec.getReferenceName();
        for (int i = 0; i < referenceSequencePrefix.length; ++i)
if (reference.startsWith(referenceSequencePrefix[i])) {
                final String geneExon = rec.getStringAttribute(GENE_NAME_TAG);
                if (geneExon != null) {
                    rec.setAttribute(GENE_NAME_TAG, genePrefix[i] + geneExon);
                    break;
                }
            }
        return rec;
    }
 

@Override
protected void processRecord(final SAMRecord rec) {
	// String recName = rec.getReadName();
	int coordinate = getCoordinate(rec);
	//Interval i = new Interval(rec.getReferenceName(), coordinate, coordinate);
	//String value = IntervalTagComparator.toString(i);
	// encode the string manually as an interval, don't use the full encoding as it's not needed!
	String value = rec.getReferenceName() + IntervalTagComparator.ENCODE_DELIMITER + coordinate;
	rec.setAttribute(this.tag, value);
	queueRecordForOutput(rec);
}
 

@Override
protected String getReferenceName(SAMRecord record) {
  return record.getReferenceName();
}
 

public void acceptRecord(final SAMRecordAndReference args) {
	mappedRecordCount++;

	final SAMRecord samRecord = args.getSamRecord();
	final ReferenceSequence referenceSequence = args.getReferenceSequence();

	final byte[] readBases = samRecord.getReadBases();
	final byte[] readQualities = samRecord.getBaseQualities();
	final byte[] refBases = referenceSequence.getBases();

	if (samRecord.getReadLength() < minReadLength) {
		smallReadCount++;
		return;
	} else if (SequenceUtil.countMismatches(samRecord, refBases, true) > Math.round(samRecord.getReadLength() * maxMismatchRate)) {
		mismatchCount++;
		return;
	}

	// We only record non-CpG C sites if there was at least one CpG in the read, keep track of
	// the values for this record and then apply to the global totals if valid
	int recordCpgs = 0;

	for (final AlignmentBlock alignmentBlock : samRecord.getAlignmentBlocks()) {
		final int blockLength = alignmentBlock.getLength();
		final int refFragmentStart = alignmentBlock.getReferenceStart() - 1;
		final int readFragmentStart = alignmentBlock.getReadStart() - 1;

		final byte[] refFragment = getFragment(refBases, refFragmentStart, blockLength);
		final byte[] readFragment = getFragment(readBases, readFragmentStart, blockLength);
		final byte[] readQualityFragment = getFragment(readQualities, readFragmentStart, blockLength);

		if (samRecord.getReadNegativeStrandFlag()) {
			// In the case of a negative strand, reverse (and complement for base arrays) the reference,
			// reads & qualities so that it can be treated as a positive strand for the rest of the process
			SequenceUtil.reverseComplement(refFragment);
			SequenceUtil.reverseComplement(readFragment);
			SequenceUtil.reverseQualities(readQualityFragment);
		}

		for (int i=0; i < blockLength-1; i++) {
			final int curRefIndex = getCurRefIndex(refFragmentStart, blockLength, i, samRecord.getReadNegativeStrandFlag());

			// Look at a 2-base window to see if we're on a CpG site, and if so check for conversion
			// (CG -> TG). We do not consider ourselves to be on a CpG site if we're on the last base of a read
			if ((SequenceUtil.basesEqual(refFragment[i], SequenceUtil.C)) &&
					(SequenceUtil.basesEqual(refFragment[i+1], SequenceUtil.G))) {
				// We want to catch the case where there's a CpG in the reference, even if it is not valid
				// to prevent the C showing up as a non-CpG C down below. Otherwise this could have been all
				// in one if statement
				if (isValidCpg(refFragment, readFragment, readQualityFragment, i)) {
					recordCpgs++;
					final CpgLocation curLocation = new CpgLocation(samRecord.getReferenceName(), curRefIndex);
					cpgTotal.increment(curLocation);
					if (SequenceUtil.isBisulfiteConverted(readFragment[i], refFragment[i])) {
						cpgConverted.increment(curLocation);
					}
				}
				i++;
			} else if (isC(refFragment[i], readFragment[i]) && isAboveCytoQcThreshold(readQualities, i) &&
					SequenceUtil.bisulfiteBasesEqual(false, readFragment[i+1], refFragment[i+1])) {
				// C base in the reference that's not associated with a CpG
				nCytoTotal++;
				if (SequenceUtil.isBisulfiteConverted(readFragment[i], refFragment[i])) {
					nCytoConverted++;
				}
			}
		}
	}

	if (recordCpgs == 0) {
		noCpgCount++;
	}
}
 

@Test
public void testMultiHit() throws IOException {
    final File unmappedSam = new File(TEST_DATA_DIR, "multihit.unmapped.sam");
    final File alignedSam = new File(TEST_DATA_DIR, "multihit.aligned.sam");
    final File merged = File.createTempFile("merged", ".sam");
    merged.deleteOnExit();

    doMergeAlignment(unmappedSam, Collections.singletonList(alignedSam),
            null, null, null, null,
            false, true, false, 1,
            "0", "1.0", "align!", "myAligner",
            true, fasta, merged,
            null, null, null, null, null, null);

    // Iterate over the merged output and gather some statistics
    final Map<String, AlignmentAccumulator> accumulatorMap = new HashMap<String, AlignmentAccumulator>();

    final SamReader reader = SamReaderFactory.makeDefault().open(merged);
    for (final SAMRecord rec : reader) {
        final String readName;
        if (!rec.getReadPairedFlag()) readName = rec.getReadName();
        else if (rec.getFirstOfPairFlag()) readName = rec.getReadName() + "/1";
        else readName = rec.getReadName() + "/2";
        AlignmentAccumulator acc = accumulatorMap.get(readName);
        if (acc == null) {
            acc = new AlignmentAccumulator();
            accumulatorMap.put(readName, acc);
        }
        if (rec.getReadUnmappedFlag()) {
            Assert.assertFalse(acc.seenUnaligned, readName);
            acc.seenUnaligned = true;
        } else {
            ++acc.numAlignments;
            if (!rec.getNotPrimaryAlignmentFlag()) {
                Assert.assertNull(acc.primaryAlignmentSequence, readName);
                acc.primaryAlignmentSequence = rec.getReferenceName();
            }
        }
    }

    // Set up expected results
    final Map<String, AlignmentAccumulator> expectedMap = new HashMap<String, AlignmentAccumulator>();
    expectedMap.put("frag_hit", new AlignmentAccumulator(false, 1, "chr1"));
    expectedMap.put("frag_multihit", new AlignmentAccumulator(false, 3, "chr2"));
    expectedMap.put("frag_no_hit", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_both_hit/1", new AlignmentAccumulator(false, 1, "chr7"));
    expectedMap.put("pair_both_hit/2", new AlignmentAccumulator(false, 1, "chr7"));
    expectedMap.put("pair_both_multihit/1", new AlignmentAccumulator(false, 3, "chr8"));
    expectedMap.put("pair_both_multihit/2", new AlignmentAccumulator(false, 3, "chr8"));
    expectedMap.put("pair_first_hit/1", new AlignmentAccumulator(false, 1, "chr1"));
    expectedMap.put("pair_first_hit/2", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_first_multihit/1", new AlignmentAccumulator(false, 3, "chr1"));
    expectedMap.put("pair_first_multihit/2", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_no_hit/1", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_no_hit/2", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_second_hit/1", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_second_hit/2", new AlignmentAccumulator(false, 1, "chr1"));
    expectedMap.put("pair_second_multihit/1", new AlignmentAccumulator(true, 0, null));
    expectedMap.put("pair_second_multihit/2", new AlignmentAccumulator(false, 3, "chr3"));

    Assert.assertEquals(accumulatorMap.size(), expectedMap.size());

    for (final Map.Entry<String, AlignmentAccumulator> entry : expectedMap.entrySet()) {
        final AlignmentAccumulator actual = accumulatorMap.get(entry.getKey());
        Assert.assertEquals(actual, entry.getValue(), entry.getKey());
    }

    assertSamValid(merged);
}
 

@Override
protected int doWork() {
	IOUtil.assertFileIsReadable(INPUT);
	IOUtil.assertFileIsWritable(OUTPUT);
	buildPatterns();

	SamReader in = SamReaderFactory.makeDefault().open(INPUT);

	SAMFileHeader fileHeader = editSequenceDictionary(in.getFileHeader().clone());
	SamHeaderUtil.addPgRecord(fileHeader, this);
	SAMFileWriter out = new SAMFileWriterFactory().makeSAMOrBAMWriter(fileHeader, true, OUTPUT);
	ProgressLogger progLog=new ProgressLogger(log);

	final boolean sequencesRemoved = fileHeader.getSequenceDictionary().getSequences().size() != in.getFileHeader().getSequenceDictionary().getSequences().size();
	
	FilteredReadsMetric m = new FilteredReadsMetric();
	
	for (final SAMRecord r : in) {
		progLog.record(r);
		if (!filterRead(r)) {
               String sequenceName = stripReferencePrefix(r.getReferenceName());
               String mateSequenceName = null;
               if (r.getMateReferenceIndex() != -1)
				mateSequenceName = stripReferencePrefix(r.getMateReferenceName());
		    if (sequencesRemoved || sequenceName != null) {
		        if (sequenceName == null)
					sequenceName = r.getReferenceName();
                   // Even if sequence name has not been edited, if sequences have been removed, need to set
                   // reference name again to invalidate reference index cache.
                   r.setReferenceName(sequenceName);
               }
               if (r.getMateReferenceIndex() != -1 && (sequencesRemoved || mateSequenceName != null)) {
                   if (mateSequenceName == null)
					mateSequenceName = r.getMateReferenceName();
                   // It's possible that the mate was mapped to a reference sequence that has been dropped from
                   // the sequence dictionary.  If so, set the mate to be unmapped.
                   if (fileHeader.getSequenceDictionary().getSequence(mateSequenceName) != null)
					r.setMateReferenceName(mateSequenceName);
				else {
                       r.setMateUnmappedFlag(true);
                       r.setMateReferenceIndex(SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX);
                       r.setMateAlignmentStart(0);
                   }
               }
			out.addAlignment(r);
			m.READS_ACCEPTED++;
		} else {
			m.READS_REJECTED++;
		}
	}
	// write the summary if the summary file is not null.
	writeSummary(this.SUMMARY, m);
       CloserUtil.close(in);
       out.close();
       FilterProgramUtils.reportAndCheckFilterResults("reads", m.READS_ACCEPTED, m.READS_REJECTED,
			PASSING_READ_THRESHOLD, log);
	return (0);
}
 

/**
 * Copied from net.sf.picard.sam.SamPairUtil. This is a more permissive
 * version of the method, which does not reset alignment start and reference
 * for unmapped reads.
 * 
 * @param rec1
 * @param rec2
 * @param header
 */
public static void setLooseMateInfo(final SAMRecord rec1, final SAMRecord rec2, final SAMFileHeader header) {
	if (rec1.getReferenceName() != SAMRecord.NO_ALIGNMENT_REFERENCE_NAME
			&& rec1.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX)
		rec1.setReferenceIndex(header.getSequenceIndex(rec1.getReferenceName()));
	if (rec2.getReferenceName() != SAMRecord.NO_ALIGNMENT_REFERENCE_NAME
			&& rec2.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX)
		rec2.setReferenceIndex(header.getSequenceIndex(rec2.getReferenceName()));

	// If neither read is unmapped just set their mate info
	if (!rec1.getReadUnmappedFlag() && !rec2.getReadUnmappedFlag()) {

		rec1.setMateReferenceIndex(rec2.getReferenceIndex());
		rec1.setMateAlignmentStart(rec2.getAlignmentStart());
		rec1.setMateNegativeStrandFlag(rec2.getReadNegativeStrandFlag());
		rec1.setMateUnmappedFlag(false);
		rec1.setAttribute(SAMTag.MQ.name(), rec2.getMappingQuality());

		rec2.setMateReferenceIndex(rec1.getReferenceIndex());
		rec2.setMateAlignmentStart(rec1.getAlignmentStart());
		rec2.setMateNegativeStrandFlag(rec1.getReadNegativeStrandFlag());
		rec2.setMateUnmappedFlag(false);
		rec2.setAttribute(SAMTag.MQ.name(), rec1.getMappingQuality());
	}
	// Else if they're both unmapped set that straight
	else if (rec1.getReadUnmappedFlag() && rec2.getReadUnmappedFlag()) {
		rec1.setMateNegativeStrandFlag(rec2.getReadNegativeStrandFlag());
		rec1.setMateUnmappedFlag(true);
		rec1.setAttribute(SAMTag.MQ.name(), null);
		rec1.setInferredInsertSize(0);

		rec2.setMateNegativeStrandFlag(rec1.getReadNegativeStrandFlag());
		rec2.setMateUnmappedFlag(true);
		rec2.setAttribute(SAMTag.MQ.name(), null);
		rec2.setInferredInsertSize(0);
	}
	// And if only one is mapped copy it's coordinate information to the
	// mate
	else {
		final SAMRecord mapped = rec1.getReadUnmappedFlag() ? rec2 : rec1;
		final SAMRecord unmapped = rec1.getReadUnmappedFlag() ? rec1 : rec2;

		mapped.setMateReferenceIndex(unmapped.getReferenceIndex());
		mapped.setMateAlignmentStart(unmapped.getAlignmentStart());
		mapped.setMateNegativeStrandFlag(unmapped.getReadNegativeStrandFlag());
		mapped.setMateUnmappedFlag(true);
		mapped.setInferredInsertSize(0);

		unmapped.setMateReferenceIndex(mapped.getReferenceIndex());
		unmapped.setMateAlignmentStart(mapped.getAlignmentStart());
		unmapped.setMateNegativeStrandFlag(mapped.getReadNegativeStrandFlag());
		unmapped.setMateUnmappedFlag(false);
		unmapped.setInferredInsertSize(0);
	}

	boolean firstIsFirst = rec1.getAlignmentStart() < rec2.getAlignmentStart();
	int insertSize = firstIsFirst ? SamPairUtil.computeInsertSize(rec1, rec2) : SamPairUtil.computeInsertSize(rec2,
			rec1);

	rec1.setInferredInsertSize(firstIsFirst ? insertSize : -insertSize);
	rec2.setInferredInsertSize(firstIsFirst ? -insertSize : insertSize);

}