function [bnet crossdat fulldat] = HMMGenSynth(T,rep)

n=3;
hs = 3; % num of hidden states or 1 if LDS

nnodes = n;% vars + 1 hidden

intra = zeros(nnodes,nnodes);
inter = zeros(nnodes,nnodes);

%NBC structure
for i = 2:nnodes
    intra(1,i)=1;
end

%AR links - trans matrices
for i = 2:nnodes
    inter(i,i)=1;
end
inter(1,1)=1;


intra
inter

ns = ones(1,nnodes*2);
ns(1)=hs; 
ns(nnodes+1)=hs;

obs = [1:nnodes];

%eclass1 = [1 2 3 4 5];
%eclass2 = [6 2 3 4 5];

ds = [1];%if HMM
%ds = [];%if LDS


bnet = mk_dbn(intra, inter, ns, 'discrete', ds, 'observed', obs);%'eclass1', eclass1, 'eclass2', eclass2, 


%for i=2:nnodes
%    bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', [0 1 1], 'cov', [.05, .1 .1], 'weights', []);
%end
%for i=nnodes+2:nnodes*2
%    bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', [0 1 1], 'cov', [.05 .1 .1], 'weights', [.5, .5 -.5]);
%end

%bnet.CPD{nnodes*2} = gaussian_CPD(bnet, i, 'mean', [0 1 1], 'cov', [.05 .1 .1], 'weights', [.5, -.5 .5]);


for i=2:nnodes
    bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', [0 1 1], 'cov', [0.3, 0.01 0.01], 'weights', []);
end
for i=nnodes+2:nnodes*2
    bnet.CPD{i} = gaussian_CPD(bnet, i, 'mean', [0 -1 1], 'cov', [0.05 0.3 0.3], 'weights', [0.8, 0.8 0.8]);
end

bnet.CPD{nnodes*2} = gaussian_CPD(bnet, i, 'mean', [0 1 1], 'cov', [0.05 0.3 0.3], 'weights', [0.8, 0.8 0.8]);

bnet.CPD{1} = tabular_CPD(bnet, 1); %if HMM
bnet.CPD{nnodes+1} = tabular_CPD(bnet, nnodes+1, [0.8 0.1 0.1 ; 0 1 0; 0 0 1]);%if HMM




%engine = smoother_engine(hmm_2TBN_inf_engine(bnet)); %rigid hmm structure
engine = smoother_engine(jtree_2TBN_inf_engine(bnet)); %no hidden variables
%engine = jtree_unrolled_dbn_inf_engine(bnet, trainlen); %slow

%[bnet2, LLtrace] = learn_params_dbn_em(engine, traindata, 'max_iter', 50, 'init_temp', 100, 'thresh', 1e-10, 'anneal', 0);


%to set the engine based on bnet2
%engine = smoother_engine(hmm_2TBN_inf_engine(bnet2));
%engine = jtree_unrolled_dbn_inf_engine(bnet, order);
%engine = smoother_engine(jtree_2TBN_inf_engine(bnet2));


fulldat=[];
crossdat=[];
for i = 1:rep
    start=cell(n,T);
    start{1,1}=1;
    newts=cell2num(sample_dbn(bnet,'evidence',start));
    fulldat=[fulldat newts];
    crossdat = [crossdat newts(:,floor(rand*T)+1)];
end

