% ---- DVH PLOT FUNCTION ----
function plot_DVH(Geometry, dose)
    % this function plots the DVHs of the given dose
    
    nfrac=1;
    
    
    colorList = jet(numel(Geometry.ROIS));
    names = {};
    
    figure
    hold on
    for roi_idx= 1:numel(Geometry.ROIS);
        % plot the histogram
        
        [dvh dosebins] = get_dvh_data(dose, Geometry.ROIS{roi_idx}.ind, nfrac);
        plot(dosebins, dvh,'Color', colorList(roi_idx, :),'LineStyle', '-','DisplayName', Geometry.ROIS{roi_idx}.name);
        names{roi_idx} = Geometry.ROIS{roi_idx}.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(['NLP optimized DVH'])
    legend(names)
    axis([0 110 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