Function getAlignmentIdentity

PopBubbles/PopBubbles.cpp:259–304

 */
template<typename It>
static float
getAlignmentIdentity(const Graph& g, vertex_descriptor t, vertex_descriptor v, It first, It last)
{
	unsigned nbranches = distance(first, last);
	vector<int> inDists(nbranches);
	transform(
	    first, last, inDists.begin(), boost::lambda::bind(getDistance, boost::cref(g), t, _1));
	vector<int> outDists(nbranches);
	transform(
	    first, last, outDists.begin(), boost::lambda::bind(getDistance, boost::cref(g), _1, v));
	vector<int> insertLens(nbranches);
	transform(
	    first,
	    last,
	    insertLens.begin(),
	    boost::lambda::bind(getDistance, boost::cref(g), t, _1) +
	        boost::lambda::bind(getLength, &g, _1) +
	        boost::lambda::bind(getDistance, boost::cref(g), _1, v));

	int max_in_overlap = -(*min_element(inDists.begin(), inDists.end()));
	assert(max_in_overlap >= 0);
	int max_out_overlap = -(*min_element(outDists.begin(), outDists.end()));
	assert(max_out_overlap >= 0);
	int min_insert_len = *min_element(insertLens.begin(), insertLens.end());
	int max_insert_len = *max_element(insertLens.begin(), insertLens.end());

	float max_identity = (float)(min_insert_len + max_in_overlap + max_out_overlap) /
	                     (max_insert_len + max_in_overlap + max_out_overlap);
	if (min_insert_len <= 0 || max_identity < opt::identity)
		return max_identity;

	vector<string> seqs(nbranches);
	transform(first, last, seqs.begin(), boost::lambda::bind(getSequence, &g, _1));
	for (unsigned i = 0; i < seqs.size(); i++) {
		// Remove the overlapping sequence.
		int n = seqs[i].size();
		int l = -inDists[i], r = -outDists[i];
		assert(n > l + r);
		seqs[i] = seqs[i].substr(l, n - l - r);
	}

	unsigned matches, consensusSize;
	tie(matches, consensusSize) = align(seqs);
	return (float)(matches + max_in_overlap + max_out_overlap) /
	       (consensusSize + max_in_overlap + max_out_overlap);
}

/** Pop the specified bubble if it is a simple bubble.
 * @return whether the bubble is popped