#include <ncbi.h>
#include <blast/blast.h>

BLAST_HSPPtr LIBCALL
evals(bcp, hp0, dblen)
	BLAST_ConfigPtr	bcp;
	BLAST_HSPPtr	hp0;
	double	dblen;
{
	register BLAST_HSPPtr	hp;
	register double	p, e;
	BLAST_KarlinBlkPtr	kbp;
	register BLAST_HSPPtr	best_hsp;

	for (hp = best_hsp = hp0; hp != NULL; hp = hp->next) {
		kbp = hp->kbp;
		e = BlastKarlinStoE(hp->score, kbp, hp->q_seg.sp->efflen, hp->s_seg.sp->efflen);
		p = BlastEtoPoissonP(hp->n, e);
		p = BlastGapDecay(p, hp->n, bcp->gapdecayrate);

		/*
		Adjust for the length of the target sequence versus
		the total length of the database
		*/
		hp->etype = VSCORE_POISSONP;
		hp->evalue = p * dblen / hp->s_seg.sp->efflen;
		hp->pvalue = -1.;

		/* Remember which HSP had the smallest expectation */
		if (hp->evalue < best_hsp->evalue)
			best_hsp = hp;
	}
	return best_hsp;
}

BLAST_HSPPtr LIBCALL
consist_evals(bcp, hp0, dblen)
	BLAST_ConfigPtr	bcp;
	BLAST_HSPPtr	hp0;
	double	dblen;
{
	register BLAST_HSPPtr	hp, best_hsp;
	register double	p, e;
	BLAST_KarlinBlkPtr	kbp;

	for (hp = best_hsp = hp0; hp != NULL; hp = hp->next) {
		kbp = hp->kbp;
		e = BlastKarlinStoE(hp->score, kbp, hp->q_seg.sp->efflen, hp->s_seg.sp->efflen);
		p = BlastEtoConsistP(hp->n, e);
		p = BlastGapDecay(p, hp->n, bcp->gapdecayrate);

		/*
		Adjust for the length of the target sequence versus
		the total length of the database
		*/
		hp->etype = VSCORE_CONSISTE;
		hp->evalue = p * dblen / hp->s_seg.sp->efflen;
		hp->pvalue = -1.;

		/* Remember which HSP had the smallest expectation */
		if (hp->evalue < best_hsp->evalue)
			best_hsp = hp;
	}
	return best_hsp;
}