<commit_before>/*
 * Copyright (C) 2009-2018 by <NAME> (<EMAIL>)
 *
 * Released under the MIT license, see LICENSE.txt
 */

#include "sonLib.h"
#include "pairwiseAlignment.h"
#include "math.h"

uint64_t getStartCoordinate(struct PairwiseAlignment *pairwiseAlignment) {
	assert(pairwiseAlignment->strand1); // This code assumes that the alignment is reported with respect
	// to the positive strand of the first sequence
	return pairwiseAlignment->start1;
}

int cmpAlignmentsFn(const void *a, const void *b) {
	const struct PairwiseAlignment *pA1 = a;
	const struct PairwiseAlignment *pA2 = b;
	return pA1->score < pA2->score ? -1 : (pA1->score > pA2->score ? 1 : 0);
}

void updateScoresToReflectMappingQualities(stList *alignments, float alpha) {
	// Create an array of the scores
	float *alignmentScores = st_calloc(stList_length(alignments), sizeof(float));
	for(uint64_t i=0; i<stList_length(alignments); i++) {
		alignmentScores[i] = ((struct PairwiseAlignment *)stList_get(alignments, i))->score;
	}

	// Calculate mapQs
	for(uint64_t i=0; i<stList_length(alignments); i++) {
		struct PairwiseAlignment *pA = stList_get(alignments, i);

		// Cut off the calculation if clearly going to be zero
		if(alpha * (alignmentScores[i] - alignmentScores[stList_length(alignments)-1]) < -10) {
			pA->score = 0.0;
		}

		else {
			// Calculate the denominator
			double z = 0.0;
			for(uint64_t j=0; j<stList_length(alignments); j++) {
				z += pow(10, alpha * (alignmentScores[j] - alignmentScores[i]));
			}
			assert(z >= 1.0);

			if(z <= 1.000001) { // Round scores to max of 60
				pA->score = 60.0;
			}
			else {
				pA->score = -10.0 * log10(1.0 - 1.0/z);
				assert(pA->score >= 0.0);
			}
		}
	}

	// Cleanup
	free(alignmentScores);
}

void reportAlignments(stList *alignments, int64_t maxAlignmentsPerSite,
		float minimumMapQValue, float alpha, FILE **fileHandleOuts) {
	// Sort by ascending score
	stList_sort(alignments, cmpAlignmentsFn);

	// Calculate the mapping qualities
	updateScoresToReflectMappingQualities(alignments, alpha);

	// Report the alignments
	for(int64_t i=0; stList_length(alignments) > 0;) {
		struct PairwiseAlignment *pairwiseAlignment = stList_pop(alignments);
		if(i < maxAlignmentsPerSite && pairwiseAlignment->score >= minimumMapQValue) {
			// Write out modified cigar
			cigarWrite(fileHandleOuts[i++], pairwiseAlignment, 0);
		}

		// Cleanup
		destructPairwiseAlignment(pairwiseAlignment);
	}
}

int main(int argc, char *argv[]) {
	/*
	 * For each alignment in the input file copy the alignment to the output file and additionally
	 * write out the alignment with the query and target sequences reversed.<commit_msg>Remove actual mapQ calculation temporarily<commit_after>/*
 * Copyright (C) 2009-2018 by <NAME> (<EMAIL>)
 *
 * Released under the MIT license, see LICENSE.txt
 */

#include "sonLib.h"
#include "pairwiseAlignment.h"
#include "math.h"

uint64_t getStartCoordinate(struct PairwiseAlignment *pairwiseAlignment) {
	assert(pairwiseAlignment->strand1); // This code assumes that the alignment is reported with respect
	// to the positive strand of the first sequence
	return pairwiseAlignment->start1;
}

int cmpAlignmentsFn(const void *a, const void *b) {
	const struct PairwiseAlignment *pA1 = a;
	const struct PairwiseAlignment *pA2 = b;
	return pA1->score < pA2->score ? -1 : (pA1->score > pA2->score ? 1 : 0);
}

void updateScoresToReflectMappingQualities(stList *alignments, float alpha) {
}

void reportAlignments(stList *alignments, int64_t maxAlignmentsPerSite,
		float minimumMapQValue, float alpha, FILE **fileHandleOuts) {
	// Sort by ascending score
	stList_sort(alignments, cmpAlignmentsFn);

	// Calculate the mapping qualities
	updateScoresToReflectMappingQualities(alignments, alpha);

	// Report the alignments
	for(int64_t i=0; stList_length(alignments) > 0;) {
		struct PairwiseAlignment *pairwiseAlignment = stList_pop(alignments);
		if(i < maxAlignmentsPerSite && pairwiseAlignment->score >= minimumMapQValue) {
			// Write out modified cigar
			cigarWrite(fileHandleOuts[i++], pairwiseAlignment, 0);
		}

		// Cleanup
		destructPairwiseAlignment(pairwiseAlignment);
	}
}

int main(int argc, char *argv[]) {
	/*
	 * For each alignment in the input file copy the alignment to the output file and additionally
	 * write out the alignment with the query and target sequences reversed.
{'fraction_non_alphanumeric': 0.11361028230433785, 'fraction_numerical': 0.02272205646086757, 'mean_word_length': 4.594351732991014, 'pattern_counts': {'":': 0, '<': 16, '<?xml version=': 0, '>': 34, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, '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': '3419695', 'n_tokens_mistral': 1575, 'n_tokens_neox': 1395, 'n_words': 394}