/************ 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); } */