function mat = mk_hmm_matrix(bnet, pdom, cdom)
% MK_HMM_MATRIX Make a matrix Pr(Y|X), where X represents all the parents, and Y all the children (assumed cond. indep.)
% mat = mk_hmm_matrix(bnet, pdom, cdom)
%
% e.g., Consider the DBN with interconnectivity i->i', j->j',k', k->i',k'
% Then transition matrix = Pr(i,j,k -> i',j',k') = Pr(i,k->i') Pr(j->j') Pr(j,k->k')


dom = [pdom cdom];
T = mk_table(dom, bnet.node_sizes(dom));
for i=1:length(cdom)
  j = cdom(i);
  e = bnet.equiv_class(j);
  % convert the CPT to a labeled table
  fam = family(bnet.adj_mat, j);
  pot = mk_table(fam, bnet.node_sizes(fam));
  pot.T = myreshape(bnet.bnodes{e}.CPT, bnet.node_sizes(fam));
  T = multiply_by_table(T, pot);
end
psize = prod(bnet.node_sizes(pdom));
csize = prod(bnet.node_sizes(cdom));
mat = reshape(T.T, [psize csize]);


