% ---- DVH PLOT FUNCTION ----
function plot_DVH(dose, optGoal)
    % this function plots the DVHs of the given dose
    
    
    
    nfrac=1;
    
    figure
    hold on
    for roi_idx=optGoal_idx
        % plot the histogram
        
        [dvh dosebins] = get_dvh_data(dose,optGoal{roi_idx}.ROI_idx,nfrac);
        plot(dosebins, dvh,'Color', optGoal{roi_idx}.dvh_col,'LineStyle', '-','DisplayName', optGoal{roi_idx}.ROI_name);
    end
    hold off
    % get % of volume above/below threshold
    num_vox = numel(optGoal{targetMinMax_idx(1)}.target);
    perc_vox_min = 100* numel(find((dose(optGoal{targetMinMax_idx(1)}.ROI_idx) - optGoal{targetMinMax_idx(1)}.target*0.95) <0)) ...
        / num_vox;
    perc_vox_max = 100* numel(find((dose(optGoal{targetMinMax_idx(2)}.ROI_idx) - optGoal{targetMinMax_idx(2)}.target*1.07) >0)) ...
        / num_vox;
    title([num2str(perc_vox_min) ' % below 0.95D_0, ' num2str(perc_vox_max) ' % above D 1.07D_0'])
    
    axis([0 100 0 100])
end
function plot_DVH_robust(dose, optGoal, optGoal_idx)
    % this function plots the DVHs of the given dose
    nfrac=1;
    
    figure
    hold on
    for goal_i=optGoal_idx
        % for all the perturbations
        for nrs_i = 1:optGoal{1}.NbrRandScenarios 
            for sss_i = 1:optGoal{1}.NbrSystSetUpScenarios % syst. setup scenarios = ss
                for rrs_i = 1:optGoal{1}.NbrRangeScenarios % range scenario = rs
                    % plot the histogram
                    ROI_idx = optGoal{goal_i}.nrs{nrs_i}.sss{sss_i}.rrs{rrs_i}.ROI_idx;
                    
                    [dvh dosebins] = get_dvh_data(dose,ROI_idx,nfrac);
                    plot(dosebins, dvh,'Color', optGoal{goal_i}.dvh_col,'LineStyle', '-','DisplayName', optGoal{goal_i}.ROI_name);
                end
            end
        end
    end
    hold off
    
end
function [dvh dosebins] = get_dvh_data(dose,roi_idx,nfrac);
    % this function calculates the data for the DVH
    dosevec = dose(roi_idx);
    
    
    [pdf dosebins] = hist(dosevec, 999);
    % clip negative bins
    pdf = pdf(dosebins >= 0);
    dosebins = dosebins(dosebins >= 0);
    dvh = fliplr(cumsum(fliplr(pdf))) / numel(dosevec) * 100;
    % append the last bin
    dosebins = [dosebins dosebins(end)+0.1];
    dvh = [dvh 0];
    
end