function h = shade_region(dir,crit,xs,dist_vals,crit_dist_vals)
% This function crosshatches the area under a distribution in the direction
% indicated by dir from the critical value in crit. Returns vector of
% handles to the crosshatched objects
% Inputs:
% dir - one value for each value in crit, with the sign of dir indicating
%       whether to shade above or below that value. For example, a dir
%       value of [-1,1] would indicate that the area below the first
%       critical value and above the second should be shaded
% crit - vector of critical values which mark the boundaries for shading
% xs - x axis values corresponding to the distribution y values
% dist_vals - y values for the distribution you want to shade regions of
%             corresponding to the x values in xs. plot(xs,dist_vals)
%             should plot the distribution you wan to shade
% crit_dist_vals - the distribution y values corresponding to the critical
%                  values listed in crit
%
% %% Mike Bindschadler 4/29/12
% %% Requires hatchfill.m

if length(dir)~=length(crit) || length(crit)~= length(crit_dist_vals)
    error('Exactly one direction must be given for each critical value and each critical value must have exactly one associated distribution value')
end

if length(xs)~=length(dist_vals)
    error('Length of x and y values don''t match for distribution')
end

%% loop over regions to shade
% If there are exactly two critical values, and we are supposed to shade
% above the greater one and below the lesser one, then only shade the
% region in between the two
if length(crit)==2 
    [maxval,i] = max(crit);
    [minval,j] = min(crit);
    if dir(i)<0 && dir(j)>=0
        % shade only in between
        mask = xs<crit(i) & xs>crit(j);
        reg_xs = [crit(j),xs(mask),crit(i),crit(i),fliplr(xs(mask)),crit(j)];
        reg_ys = [crit_dist_vals(j),dist_vals(mask),crit_dist_vals(i),0,zeros(size(xs(mask))),0];
        ph = patch(reg_xs,reg_ys,'r');
        h = hatchfill(ph,'cross'); 
        set(h,'Color','red');      
        return
    end
end

for i = 1:length(dir)
    if dir(i)<0
        mask = xs<crit(i); % Shade below value
        reg_xs = [xs(mask),crit(i),crit(i),fliplr(xs(mask))];
        reg_ys = [dist_vals(mask),crit_dist_vals(i),0,zeros(size(xs(mask)))];
    else
        mask = xs>crit(i); % Shade above value
        reg_xs = [crit(i),xs(mask),fliplr(xs(mask)),crit(i)];
        reg_ys = [crit_dist_vals(i),dist_vals(mask),zeros(size(xs(mask))),0];
    end
    ph = patch(reg_xs,reg_ys,'r');
    h(i) = hatchfill(ph,'cross'); set(h(i),'Color','red');
end