<commit_before>/**
 *   Copyright 2013 <NAME>
 *
 *   Licensed under the Apache License, Version 2.0 (the "License");
 *   you may not use this file except in compliance with the License.
 *   You may obtain a copy of the License at
 *
 *       http://www.apache.org/licenses/LICENSE-2.0
 *
 *   Unless required by applicable law or agreed to in writing, software
 *   distributed under the License is distributed on an "AS IS" BASIS,
 *   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *   See the License for the specific language governing permissions and
 *   limitations under the License.
 */

package com.roche.heatseq.process;

import java.io.PrintWriter;
import java.util.ArrayList;
import java.util.List;

import net.sf.samtools.SAMFileHeader;
import net.sf.samtools.SAMRecord;

import org.slf4j.Logger;
import org.slf4j.LoggerFactory;

import com.roche.heatseq.objects.IReadPair;
import com.roche.heatseq.objects.Probe;
import com.roche.heatseq.objects.ReadPair;
import com.roche.mapping.SAMRecordUtil;
import com.roche.sequencing.bioinformatics.common.alignment.CigarString;
import com.roche.sequencing.bioinformatics.common.alignment.IAlignmentScorer;
import com.roche.sequencing.bioinformatics.common.alignment.NeedlemanWunschGlobalAlignment;
import com.roche.sequencing.bioinformatics.common.sequence.ICode;
import com.roche.sequencing.bioinformatics.common.sequence.ISequence;
import com.roche.sequencing.bioinformatics.common.sequence.IupacNucleotideCode;
import com.roche.sequencing.bioinformatics.common.sequence.IupacNucleotideCodeSequence;
import com.roche.sequencing.bioinformatics.common.sequence.Strand;
import com.roche.sequencing.bioinformatics.common.utils.StringUtil;

/**
 * Extends reads for a probe to the primer locations
 * 
 */
public final class ExtendReadsToPrimer {
	private static final Logger logger = LoggerFactory.getLogger(ExtendReadsToPrimer.class);

	private static final int PRIMER_ALIGNMENT_BUFFER = 10;
	private static final int PRIMER_ACCEPTANCE_BUFFER = 2;

	/**
	 * We just use static methods from this class
	 */
	private ExtendReadsToPrimer() {
		throw new AssertionError();
	}

	/**
	 * Extends reads for a probe to the primer locations
	 * 
	 * @param readPair
	 * @return UniqueProbeRepresentativeData with extended reads or null if it could not be extended
	 */
	private static IReadPair extendReadPair(Probe probe, IReadPair readPair, IAlignmentScorer alignmentScorer) {
		return extendReadPair(readPair.getUid(), probe, readPair.getSamHeader(), readPair.getSequenceName(), readPair.getReadName(), readPair.getReadGroup(),
				new IupacNucleotideCodeSequence(readPair.getSequenceOne()), readPair.getSequenceOneQualityString(), new IupacNucleotideCodeSequence(readPair.getSequenceTwo()),
				readPair.getSequenceTwoQualityString(), readPair.getOneMappingQuality(), readPair.getTwoMappingQuality(), alignmentScorer);
	}

	/**
	 * 
	 * @param uid
	 * @param probe
	 * @param samHeader
	 * @param sequenceName
	 * @param readName
	 * @param readGroup
	 * @param sequenceOne
	 * @param sequenceOneQualityString
	 * @param sequenceTwo
	 * @param sequenceTwoQualityString
	 * @param oneMappingQuality
	 * @param twoMappingQuality
	 * @return readPair that has been extended to the primers (the primers are not included in the new alignment)
	 */
	public static IReadPair extendReadPair(String uid, Probe probe, SAMFileHeader samHeader, String sequenceName, String readName, String readGroup, ISequence sequenceOne,
			String sequenceOneQualityString, ISequence sequenceTwo, String sequenceTwoQualityString, int oneMappingQuality, int twoMappingQuality, IAlignmentScorer alignmentScorer) {
		IReadPair extendedReadPair = null;

		try {
			sequenceTwo = sequenceTwo.getCompliment();

			ISequence captureTargetSequence = probe.getCaptureTargetSequence();
			ISequence extensionPrimer = probe.getExtensionPrimerSequence();
			ISequence ligationPrimer = probe.getLigationPrimerSequence();
			boolean readOneIsOnReverseStrand = probe.getProbeStrand() == Strand.REVERSE;
			boolean readTwoIsOnReverseStrand = !readOneIsOnReverseStrand;

			int primerReferencePositionAdjacentToSequence = probe.getExtensionPrimerStop();
			if (readOneIsOnReverseStrand) {
				primerReferencePositionAdjacentToSequence = probe.getExtensionPrimerStart();
			}
			ReadExtensionDetails readOneExtensionDetails = calculateDetailsForReadExtensionToPrimer(extensionPrimer, primerReferencePositionAdjacentToSequence, captureTargetSequence, sequenceOne,
					false, readOneIsOnReverseStrand, alignmentScorer);

			if (readOneExtensionDetails != null) {
				String readOneExtendedSequence = sequenceOne.subSequence(readOneExtensionDetails.getReadStart(), sequenceOne.size()).toString();
				String readOneExtendedBaseQualities = sequenceOneQualityString.substring(readOneExtensionDetails.getReadStart(), sequenceOneQualityString.length()).toString();
				int readOneReferenceLength = probe.getCaptureTargetSequence().size();
				primerReferencePositionAdjacentToSequence = probe.getLigationPrimerStart();
				if (readOneIsOnReverseStrand) {
					readOneExtendedSequence = StringUtil.reverse(readOneExtendedSequence);
					readOneExtendedBaseQualities = StringUtil.reverse(readOneExtendedBaseQualities);
					readOneReferenceLength = -readOneReferenceLength;
					primerReferencePositionAdjacentToSequence = probe.getLigationPrimerStop();
				}

				SAMRecord readOneExtendedRecord = extendRecord(samHeader, readName, readGroup, readOneIsOnReverseStrand, readOneExtensionDetails.getAlignmentCigarString(),
						readOneExtensionDetails.getMismatchDetailsString(), readOneExtensionDetails.getAlignmentStartInReference(), readOneExtendedSequence, readOneExtendedBaseQualities,
						sequenceName, oneMappingQuality, probe.getCaptureTargetSequence().size(), uid);

				ReadExtensionDetails readTwoExtensionDetails = calculateDetailsForReadExtensionToPrimer(ligationPrimer, primerReferencePositionAdjacentToSequence, captureTargetSequence, sequenceTwo,
						true, readTwoIsOnReverseStrand, alignmentScorer);

				if (readTwoExtensionDetails != null) {
					String readTwoExtendedBaseQualities = sequenceTwoQualityString.substring(readTwoExtensionDetails.getReadStart(), sequenceTwoQualityString.length()).toString();
					String readTwoExtendedSequence = sequenceTwo.subSequence(readTwoExtensionDetails.getReadStart(), sequenceTwo.size()).toString();
					int readTwoReferenceLength = probe.getCaptureTargetSequence().size();
					if (readTwoIsOnReverseStrand) {
						readTwoExtendedSequence = StringUtil.reverse(readTwoExtendedSequence);
						readTwoExtendedBaseQualities = StringUtil.reverse(readTwoExtendedBaseQualities);
						readTwoReferenceLength = -readTwoReferenceLength;
					}
					SAMRecord readTwoExtendedRecord = extendRecord(samHeader, readName, readGroup, readTwoIsOnReverseStrand, readTwoExtensionDetails.getAlignmentCigarString(),
							readTwoExtensionDetails.getMismatchDetailsString(), readTwoExtensionDetails.getAlignmentStartInReference(), readTwoExtendedSequence, readTwoExtendedBaseQualities,
							sequenceName, twoMappingQuality, readTwoReferenceLength, uid);

					SAMRecordUtil.setSAMRecordsAsPair(readOneExtendedRecord, readTwoExtendedRecord);

					extendedReadPair = new ReadPair(readOneExtendedRecord, readTwoExtendedRecord, uid);

				}
			}
		} catch (Exception e) {
			logger.warn(e.getMessage(), e);
		}

		return extendedReadPair;
	}

	private static ReadExtensionDetails calculateDetailsForReadExtensionToPrimer(ISequence primerSequence, int primerReferencePositionAdjacentToSequence, ISequence captureSequence,
			ISequence readSequence, boolean isLigationPrimer, boolean isReversed, IAlignmentScorer alignmentScorer) {
		ReadExtensionDetails readExtensionDetails = null;

		if (isLigationPrimer) {
			primerSequence = primerSequence.getReverse();
			captureSequence = captureSequence.getReverse();
		}

		Integer primerEndIndexInRead = getPrimerEndIndexInRead(primerSequence, readSequence, alignmentScorer);
		boolean primerAlignedSuccesfully = (primerEndIndexInRead != null) && (primerEndIndexInRead >= 0) && (primerEndIndexInRead < readSequence.size());

		if (primerAlignedSuccesfully) {
			int captureTargetStartIndexInRead = primerEndIndexInRead + 1;
			ISequence readWithoutPrimer = readSequence.subSequence(captureTargetStartIndexInRead, readSequence.size());
			NeedlemanWunschGlobalAlignment readAlignmentWithReference = new NeedlemanWunschGlobalAlignment(captureSequence, readWithoutPrimer, alignmentScorer);
			boolean readAlignedSuccesfully = readAlignmentWithReference.getLengthNormalizedAlignmentScore() > 0;

			if (readAlignedSuccesfully) {
				CigarString cigarString = readAlignmentWithReference.getCigarString();
				String mismatchDetailsString = readAlignmentWithReference.getMismatchDetailsString();
				int alignmentStartInReference = 0;

				// TODO <NAME> 7/23/2013 need to examine if this logic is correct
				if (isReversed) {
					// ISequence readAlignment = readAlignmentWithReference.getAlignmentPair().getQueryAlignment();
					ISequence referenceAlignment = readAlignmentWithReference.getAlignmentPair().getReferenceAlignment();
					cigarString = readAlignmentWithReference.getReverseCigarString();
					mismatchDetailsString = readAlignmentWithReference.getReverseMismatchDetailsString();
					alignmentStartInReference = primerReferencePositionAdjacentToSequence - referenceAlignment.size();
				} else {
					alignmentStartInReference = primerReferencePositionAdjacentToSequence;
				}<commit_msg>fixing issue with location and sequence of results in bam<commit_after>/**
 *   Copyright 2013 <NAME>
 *
 *   Licensed under the Apache License, Version 2.0 (the "License");
 *   you may not use this file except in compliance with the License.
 *   You may obtain a copy of the License at
 *
 *       http://www.apache.org/licenses/LICENSE-2.0
 *
 *   Unless required by applicable law or agreed to in writing, software
 *   distributed under the License is distributed on an "AS IS" BASIS,
 *   WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *   See the License for the specific language governing permissions and
 *   limitations under the License.
 */

package com.roche.heatseq.process;

import java.io.PrintWriter;
import java.util.ArrayList;
import java.util.List;

import net.sf.samtools.SAMFileHeader;
import net.sf.samtools.SAMRecord;

import org.slf4j.Logger;
import org.slf4j.LoggerFactory;

import com.roche.heatseq.objects.IReadPair;
import com.roche.heatseq.objects.Probe;
import com.roche.heatseq.objects.ReadPair;
import com.roche.mapping.SAMRecordUtil;
import com.roche.sequencing.bioinformatics.common.alignment.CigarString;
import com.roche.sequencing.bioinformatics.common.alignment.IAlignmentScorer;
import com.roche.sequencing.bioinformatics.common.alignment.NeedlemanWunschGlobalAlignment;
import com.roche.sequencing.bioinformatics.common.sequence.ICode;
import com.roche.sequencing.bioinformatics.common.sequence.ISequence;
import com.roche.sequencing.bioinformatics.common.sequence.IupacNucleotideCode;
import com.roche.sequencing.bioinformatics.common.sequence.IupacNucleotideCodeSequence;
import com.roche.sequencing.bioinformatics.common.sequence.Strand;
import com.roche.sequencing.bioinformatics.common.utils.StringUtil;

/**
 * Extends reads for a probe to the primer locations
 * 
 */
public final class ExtendReadsToPrimer {
	private static final Logger logger = LoggerFactory.getLogger(ExtendReadsToPrimer.class);

	private static final int PRIMER_ALIGNMENT_BUFFER = 10;
	private static final int PRIMER_ACCEPTANCE_BUFFER = 2;

	/**
	 * We just use static methods from this class
	 */
	private ExtendReadsToPrimer() {
		throw new AssertionError();
	}

	/**
	 * Extends reads for a probe to the primer locations
	 * 
	 * @param readPair
	 * @return UniqueProbeRepresentativeData with extended reads or null if it could not be extended
	 */
	private static IReadPair extendReadPair(Probe probe, IReadPair readPair, IAlignmentScorer alignmentScorer) {
		return extendReadPair(readPair.getUid(), probe, readPair.getSamHeader(), readPair.getSequenceName(), readPair.getReadName(), readPair.getReadGroup(),
				new IupacNucleotideCodeSequence(readPair.getSequenceOne()), readPair.getSequenceOneQualityString(), new IupacNucleotideCodeSequence(readPair.getSequenceTwo()),
				readPair.getSequenceTwoQualityString(), readPair.getOneMappingQuality(), readPair.getTwoMappingQuality(), alignmentScorer);
	}

	/**
	 * 
	 * @param uid
	 * @param probe
	 * @param samHeader
	 * @param sequenceName
	 * @param readName
	 * @param readGroup
	 * @param sequenceOne
	 * @param sequenceOneQualityString
	 * @param sequenceTwo
	 * @param sequenceTwoQualityString
	 * @param oneMappingQuality
	 * @param twoMappingQuality
	 * @return readPair that has been extended to the primers (the primers are not included in the new alignment)
	 */
	public static IReadPair extendReadPair(String uid, Probe probe, SAMFileHeader samHeader, String sequenceName, String readName, String readGroup, ISequence sequenceOne,
			String sequenceOneQualityString, ISequence sequenceTwo, String sequenceTwoQualityString, int oneMappingQuality, int twoMappingQuality, IAlignmentScorer alignmentScorer) {
		IReadPair extendedReadPair = null;

		try {
			sequenceTwo = sequenceTwo.getCompliment();

			ISequence captureTargetSequence = probe.getCaptureTargetSequence();
			ISequence extensionPrimer = probe.getExtensionPrimerSequence();
			ISequence ligationPrimer = probe.getLigationPrimerSequence();
			boolean readOneIsOnReverseStrand = probe.getProbeStrand() == Strand.REVERSE;
			boolean readTwoIsOnReverseStrand = !readOneIsOnReverseStrand;

			int primerReferencePositionAdjacentToSequence = probe.getExtensionPrimerStop();
			if (readOneIsOnReverseStrand) {
				primerReferencePositionAdjacentToSequence = probe.getExtensionPrimerStart();
			}
			ReadExtensionDetails readOneExtensionDetails = calculateDetailsForReadExtensionToPrimer(extensionPrimer, primerReferencePositionAdjacentToSequence, captureTargetSequence, sequenceOne,
					false, readOneIsOnReverseStrand, alignmentScorer);

			if (readOneExtensionDetails != null) {
				ISequence readOneExtendedSequence = sequenceOne.subSequence(readOneExtensionDetails.getReadStart(), sequenceOne.size());
				String readOneExtendedBaseQualities = sequenceOneQualityString.substring(readOneExtensionDetails.getReadStart(), sequenceOneQualityString.length()).toString();
				int readOneReferenceLength = probe.getCaptureTargetSequence().size();
				primerReferencePositionAdjacentToSequence = probe.getLigationPrimerStart();
				if (readOneIsOnReverseStrand) {
					readOneExtendedSequence = readOneExtendedSequence.getReverseCompliment();
					readOneExtendedBaseQualities = StringUtil.reverse(readOneExtendedBaseQualities);
					readOneReferenceLength = -readOneReferenceLength;
					primerReferencePositionAdjacentToSequence = probe.getLigationPrimerStop();
				}

				SAMRecord readOneExtendedRecord = extendRecord(samHeader, readName, readGroup, readOneIsOnReverseStrand, readOneExtensionDetails.getAlignmentCigarString(),
						readOneExtensionDetails.getMismatchDetailsString(), readOneExtensionDetails.getAlignmentStartInReference(), readOneExtendedSequence.toString(), readOneExtendedBaseQualities,
						sequenceName, oneMappingQuality, probe.getCaptureTargetSequence().size(), uid);

				ReadExtensionDetails readTwoExtensionDetails = calculateDetailsForReadExtensionToPrimer(ligationPrimer, primerReferencePositionAdjacentToSequence, captureTargetSequence, sequenceTwo,
						true, readTwoIsOnReverseStrand, alignmentScorer);

				if (readTwoExtensionDetails != null) {
					String readTwoExtendedBaseQualities = sequenceTwoQualityString.substring(readTwoExtensionDetails.getReadStart(), sequenceTwoQualityString.length());
					ISequence readTwoExtendedSequence = sequenceTwo.subSequence(readTwoExtensionDetails.getReadStart(), sequenceTwo.size()).getCompliment();
					int readTwoReferenceLength = probe.getCaptureTargetSequence().size();
					if (readTwoIsOnReverseStrand) {
						readTwoExtendedSequence = readTwoExtendedSequence.getReverseCompliment();
						readTwoExtendedBaseQualities = StringUtil.reverse(readTwoExtendedBaseQualities);
						readTwoReferenceLength = -readTwoReferenceLength;
					}
					SAMRecord readTwoExtendedRecord = extendRecord(samHeader, readName, readGroup, readTwoIsOnReverseStrand, readTwoExtensionDetails.getAlignmentCigarString(),
							readTwoExtensionDetails.getMismatchDetailsString(), readTwoExtensionDetails.getAlignmentStartInReference(), readTwoExtendedSequence.toString(),
							readTwoExtendedBaseQualities, sequenceName, twoMappingQuality, readTwoReferenceLength, uid);

					SAMRecordUtil.setSAMRecordsAsPair(readOneExtendedRecord, readTwoExtendedRecord);

					extendedReadPair = new ReadPair(readOneExtendedRecord, readTwoExtendedRecord, uid);

				}
			}
		} catch (Exception e) {
			logger.warn(e.getMessage(), e);
		}

		return extendedReadPair;
	}

	private static ReadExtensionDetails calculateDetailsForReadExtensionToPrimer(ISequence primerSequence, int primerReferencePositionAdjacentToSequence, ISequence captureSequence,
			ISequence readSequence, boolean isLigationPrimer, boolean isReversed, IAlignmentScorer alignmentScorer) {
		ReadExtensionDetails readExtensionDetails = null;

		if (isLigationPrimer) {
			primerSequence = primerSequence.getReverse();
			captureSequence = captureSequence.getReverse();
		}

		Integer primerEndIndexInRead = getPrimerEndIndexInRead(primerSequence, readSequence, alignmentScorer);
		boolean primerAlignedSuccesfully = (primerEndIndexInRead != null) && (primerEndIndexInRead >= 0) && (primerEndIndexInRead < readSequence.size());

		if (primerAlignedSuccesfully) {
			int captureTargetStartIndexInRead = primerEndIndexInRead + 1;
			ISequence readWithoutPrimer = readSequence.subSequence(captureTargetStartIndexInRead, readSequence.size());
			NeedlemanWunschGlobalAlignment readAlignmentWithReference = new NeedlemanWunschGlobalAlignment(captureSequence, readWithoutPrimer, alignmentScorer);
			boolean readAlignedSuccesfully = readAlignmentWithReference.getLengthNormalizedAlignmentScore() > 0;

			if (readAlignedSuccesfully) {
				CigarString cigarString = readAlignmentWithReference.getCigarString();
				String mismatchDetailsString = readAlignmentWithReference.getMismatchDetailsString();
				int alignmentStartInReference = 0;

				// TODO <NAME> 7/23/2013 need to examine if this logic is correct
				if (isReversed) {
					ISequence referenceAlignment = readAlignmentWithReference.getAlignmentPair().getReferenceAlignmentWithoutTerminalInserts();
					cigarString = readAlignmentWithReference.getReverseCigarString();
					mismatchDetailsString = readAlignmentWithReference.getReverseMismatchDetailsString();
					alignmentStartInReference = primerReferencePositionAdjacentToSequence - referenceAlignment.size();
				} else {
					alignmentStartInReference = primerReferencePositionAdjacentToSequence + 1;
				}
{'fraction_non_alphanumeric': 0.07606181665264929, 'fraction_numerical': 0.00257569386038688, 'mean_word_length': 6.339891975308642, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2657360', 'n_tokens_mistral': 5819, 'n_tokens_neox': 4999, 'n_words': 1059}