Java SAMRecord.setReadString方法代码示例

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
private SAMRecord createSAMRecord(String read, String cigar, String md, int NM, boolean reverse) {
	SAMFileHeader header = new SAMFileHeader();
	List<SAMProgramRecord> PRs = new ArrayList<>();
	PRs.add(new SAMProgramRecord("Bowtie"));
	header.setProgramRecords(PRs);

	SAMRecord record = new SAMRecord(header);
	record.setReadString(read);
	record.setAttribute("NM", NM); //number of mismatches
	record.setCigarString(cigar);
	record.setAttribute("MD", md);

	if (reverse) {
		record.setReadNegativeStrandFlag(true);
	}

	return record;
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
public static SAMRecord convertToPairedEndFragment(SAMRecord rec) {
	int insertSize = Math.abs(rec.getInferredInsertSize());
	//TODO: possible to remove (insertSize > 0) clause? What does an insert size of 0 mean?
	if (rec.getReadPairedFlag() && rec.getProperPairFlag() && !rec.getReadUnmappedFlag() 
			&& (insertSize <= MAX_INSERT) && (rec.getAlignmentStart() < rec.getMateAlignmentStart())
			&& rec.getReferenceName() != "chrM") // current paired end representation doesn't do well with circular chromosomes
	{
		//We need to get the full fragment contained by read.getStart to pair.getEnd
		int readEnd=rec.getAlignmentEnd();
		int mateStart=rec.getMateAlignmentStart();
					
		//int extension = insertSize - rec.getReadLength();
		int extension=(mateStart-readEnd)+rec.getReadLength();
		if (extension <= MAX_INSERT) {
			String newRead = rec.getReadString() + StringUtils.repeat("N", extension);
			rec.setReadString(newRead);
			String newQual = StringUtils.repeat("A", newRead.length());
			if(!rec.getBaseQualityString().equals("*")) newQual = rec.getBaseQualityString() + StringUtils.repeat("A", extension);
			rec.setBaseQualityString(newQual);
			String newCigar=newRead.length()+"M";
			rec.setCigarString(newCigar);

			// Change attributes to represent single read
			rec.setMateReferenceName("*");
			rec.setMateAlignmentStart(0);
			rec.setFirstOfPairFlag(false);
			rec.setMateNegativeStrandFlag(false);
			rec.setMateUnmappedFlag(false);
			rec.setProperPairFlag(false);
			rec.setReadPairedFlag(false);
			rec.setSecondOfPairFlag(false);
		}
	} else {
		rec = null;
	}
	return rec;
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
private boolean trim(SAMRecord record, int trimMismatches) {
	record.setAttribute("XT", "true");
	String sequenceCT = record.getReadString();

	int mismatches = (Integer) record.getAttribute("NM");

	if (mismatches >= trimMismatches) {

		boolean isReverse = record.getReadNegativeStrandFlag();


		int currentMismatches = 0;
		int trimPoint = isReverse ? sequenceCT.length() : 0;

		CigarIterator cigarIterator =
				new CigarIterator(record.getCigarString(), isReverse);

		MDTagIterator mdIterator = new MDTagIterator((String) record.getAttribute("MD"), isReverse);

		CigarOperator lastCigarOperator = null;
		while (cigarIterator.hasNext() && currentMismatches < trimMismatches) {
			CigarOperator operator = cigarIterator.next();

			switch (operator) {
				case M:
					MDTagOperator mdOperator = mdIterator.next();
					if (mdOperator == MDTagOperator.VARIANT) {
						currentMismatches++;
						if (currentMismatches < trimMismatches) {
							trimPoint = trimPoint + (isReverse ? -1 : 1);
						}
					} else if (mdOperator == MDTagOperator.SEQUENCE_MATCH) {
						// do nothing
						trimPoint = trimPoint + (isReverse ? -1 : 1);
					} else {
						//should not happen
						return false; // do not trim
					}

					break;
				case I:
					currentMismatches++;
					if (currentMismatches < trimMismatches) {
						trimPoint = trimPoint + (isReverse ? -1 : 1);
					}
					break;
				case D:
					currentMismatches++;
					mdIterator.next();
					break;
				default:
					return false; //do not trim this read, since it has non-supported CIGAR operations
			}
			lastCigarOperator = operator;
		}

		if (currentMismatches >= trimMismatches) {
			if (isReverse) {
				//trimPoint ++; //trim the last mismatch
				record.setReadString(record.getReadString().substring(0, trimPoint).replaceAll(".", "" +
						TRIMMED_BASE)
						+ record.getReadString().substring(trimPoint));
			} else {
				//trimPoint --; //trim the last mismatch
				record.setReadString(record.getReadString().substring(0, trimPoint) + record.getReadString()
						.substring
								(trimPoint).replaceAll(".", "" + TRIMMED_BASE));
			}
			return true;
		}

	}
	return false;

}