% 
% Input:	
%		haplocells -  haplotype_string array in cells as cgtcangcttac
%		pids - pattern id (not unique) 501, 502, 502, ...
% Output:	
%		ctables - contingency tables snplength cells of #patterns x 2		
%
%		each ctable 
%			a	g	
%		500	4	6
%		501	5	5
%		502	3	2	with missing
%
%		a-97, c-99, g-103, n-110, t-116
% Modified:
% Date: 5/5/2006


function ctables = ctablegen(haplocells, pids)
	
	upatns = unique(pids);	%find unique patterns 501,502,... for ctable rows
	ulen = length(upatns);

	% create numeric matrix of raw snp (acgt...)

	sm = length(haplocells);	%number of snps in raw data
	sn = length(haplocells{1});	%length of each snp
	snpmat = zeros(sm, sn);

	for i=1:sm
		for j=1:sn
			snpmat(i,j) = haplocells{i}(j);
		end
	end

	ctables = cell(sn, 1);

	for i=1:sn
		ctables{i} = zeros(ulen+1, 2);
%		ctables{i} = zeros(ulen, 2);
	end

	for i=1:sn
		ctable = zeros(ulen+1, 2);
		% find snp base (agn or ctn) for each column-snp
		snpcol = snpmat(:,i);
		uacgt = unique(snpcol);

		c1 = uacgt(1);
		c2 = uacgt(2);
		if find(uacgt==110) == 2
			c2 = uacgt(3);
		end

		ctable(1,1) = c1;
		ctable(1,2) = c2;

		for k=1:ulen
			upatn = upatns(k);			%each unique pattern
			snpupatn = snpcol(find(pids==upatn));	%corresponding snp values
			ctable(k+1, 1) = length(find(snpupatn==c1));
			ctable(k+1, 2) = length(find(snpupatn==c2));
		end
		
		ctable;
		ctables{i} = ctable;
	end