/************ make_gametes.c  *******************************************
*
*
*****************************************************************************/
#include "my.h"

#define STATE1 '0'
#define STATE2 '1'

struct node{
	int abv;
	int ndes;
	float time;
/*	unsigned des[4];  */
	};

	int
make_gametes(nsam, ptree,tt,newsites, ns, list )
	int nsam, newsites, ns;
	struct node *ptree;
	char **list;
	double tt ;
{
	int  tip, j,  node ;

	for(  j=ns; j< ns+newsites ;  j++ ) {
		node = pickb(  nsam, ptree, tt);
		for( tip=0; tip < nsam ; tip++) {
		   if( tdesn(ptree, tip, node) ) list[tip][j] = STATE1 ;
		   else list[tip][j] = STATE2 ;
		   }
		}
}


/***  ttime.c : Returns the total time in the tree, *ptree, with nsam tips. **/

	double
ttime( ptree, nsam)
	struct node *ptree;
	int nsam;
{
	double t;
	int i;

	t = (ptree + 2*nsam-2) -> time ;
	for( i=nsam; i< 2*nsam-1 ; i++)
		t += (ptree + i)-> time ;
	return(t);
}

/***  pickb : returns a random branch from the tree. The probability of picking
              a particular branch is proportional to its duration. tt is total
	      time in tree.   ****/

	int
pickb(nsam, ptree, tt)
	int nsam;
	struct node *ptree;
	double tt;
{
	double x, y, ran1(), dummy;
	int i;

	dummy = ran1();
	x = dummy * tt;

	for( i=0, y=0; i < 2*nsam-2 ; i++) {
		y += (ptree + (ptree+i)->abv )->time - (ptree+i)->time ;
		if( y >= x ) return( i ) ;
		}
	return( i );
}

/****  tdesn : returns 1 if tip is a descendant of node in *ptree, otherwise 0. **/

	int
tdesn(ptree, tip, node )
	struct node *ptree;
	int tip, node;
{
	int k;

	for( k= tip ; k < node ; k = (ptree+k)->abv ) ;
	if( k==node ) return(1);
	else return(0);
}

/*** poisso : returns a poisson distributed random variable with mean u. **/
/*
	int
poisso(u)
	double u;
{
	double  cump, ru, ran1(), p;
	int i=1;

	ru = ran1();
	p = exp(-u);
	if( ru < p) return(0);
	cump = p;
	
	while( ru > ( cump += (p *= u/i ) ) )
		i++;
	return(i);
}
*/

