function query_prob = BR_KB_query_prob(KB,vtable,x,query)
% BR_KB_query_prob - determine probability of KB query
% On input:
%     KB (KB struct): knowledge base CNF
%     vtable (kx7 array): variable table
%     x (kx1 vector): variable probabilities
%     query (1xp vector): p is 1, 2 or 3 (a disjunction)
% On output:
%     query_prob (float): probability of query
% Call:
%     p = BR_KB_query_prob(KB1,vtable1,x1,[1,2,-3]);
% Author:
%     T. Henderson
%     UU
%     Summer 2018
%

query_prob = 0;
F = [];

if isempty(query)
    return
end

if length(query)==1
    atom1 = abs(query);
    if query>0
        query_prob = x(atom1);
    else
        query_prob = 1 - x(atom1);
    end
end

if length(query)==2
    L1 = query(1);
    L2 = query(2);
    atom1 = abs(L1);
    atom2 = abs(L2);
    if L1>0
        p1 = x(atom1);
    else
        p1 = 1 - x(atom1);
    end
    if L2>0
        p2 = x(atom2);
    else
        p2 = 1 - x(atom2);
    end
    index = BR_get_vtable(vtable,L1,L2,0,0,BR_parity(L1,L2,0),0);
    if ~isempty(index)
        query_prob = p1 + p2 - x(index)*p2;
    else
        F = ['F(1) = ',num2str(p1),' + ',num2str(p2),' - x(2)*',...
            num2str(p2),' - x(1);'];
    end
end

if ~isempty(F)&length(query)==2
    fd = fopen('tjqwa.m','w');
    fprintf(fd,'function F = %s(x)\n','tjqwa');
    fprintf(fd,'%%\n\n');
    fprintf(fd,'%s\n',F);
    fclose(fd);
    [xg,eg,xt,et] = BR_KB_Jacobian(KB,'tjqwa',0.0001,0.0001,20000,...
        0.5*ones(2,1));
    query_prob = xg(1);
end

if length(query)==3
    L1 = query(1);
    L2 = query(2);
    L3 = query(3);
    atom1 = abs(L1);
    atom2 = abs(L2);
    atom3 = abs(L3);
    if L1>0
        p1 = x(atom1);
    else
        p1 = 1 - x(atom1);
    end
    if L2>0
        p2 = x(atom2);
    else
        p2 = 1 - x(atom2);
    end
    if L3>0
        p3 = x(atom3);
    else
        p3 = 1 - x(atom3);
    end
    num_unknowns = 1;
    index12 = BR_get_vtable(vtable,L1,L2,0,0,BR_parity(L1,L2,0),0);
    index13 = BR_get_vtable(vtable,L1,L3,0,0,BR_parity(L1,L3,0),0);
    index23 = BR_get_vtable(vtable,L2,L3,0,0,BR_parity(L2,L3,0),0);
    index123 = BR_get_vtable(vtable,L1,L2,L3,0,BR_parity(L1,L2,L3),1);
    if ~isempty(index12)
        s12 = num2str(x(index12));
    else
        num_unknowns = num_unknowns + 1;
        s12 = ['x(',num2str(num_unknowns),')'];
    end
    if ~isempty(index13)
        s13 = num2str(x(index13));
    else
        num_unknowns = num_unknowns + 1;
        s13 = ['x(',num2str(num_unknowns),')'];
    end
    if ~isempty(index23)
        s23 = num2str(x(index23));
    else
        num_unknowns = num_unknowns + 1;
        s23 = ['x(',num2str(num_unknowns),')'];
    end
    if ~isempty(index123)
        s123 = num2str(x(index123));
    else
        num_unknowns = num_unknowns + 1;
        s123 = ['x(',num2str(num_unknowns),')'];
    end
    F = ['F(1) = -x(1) + ',num2str(p1),' + ',num2str(p2),' + ',...
        num2str(p3),' - ',s12,'*',num2str(p2),' - ',s13,'*',num2str(p3),...
        ' - ',s23,'*',num2str(p3),' + ',s123,'*',s23,'*',num2str(p3),';'];
end

if ~isempty(F)&length(query)==3
    fd = fopen('tjqwa.m','w');
    fprintf(fd,'function F = %s(x)\n','tjqwa');
    fprintf(fd,'%%\n\n');
    fprintf(fd,'%s\n',F);
    fclose(fd);
    [xg,eg,xt,et] = BR_KB_Jacobian(KB,'tjqwa',0.0001,0.0001,20000,...
        0.5*ones(num_unknowns,1));
    query_prob = xg(1);
end

tch = 0;
