function mod_UQLab_plotSeriesPred(Sample, Data, plotType, varargin)
% function plotSeriesPred from UQLab File
%         uq_display_uq_inversion.m (v.1.3) 
% with some modifications  by O. Klein, WIAS, Berlin
%% ------------------------------------------------------------------------
%PLOTSERIESPRED creates a plot of serial (non-scalar) predictions.

%% Parse inputs
%% O. Klein: Modified to get values for StartInd and EndInd
%% Start of block with modifications
 StartInd = 1;
    EndInd = size(Data.y,2);
if nargin == 3
    % Point estimate and Index boundaries not available
    pointEstimatePredFlag = false;
  
else
    % Point estimate available
    pointEstimatePredFlag = true;
    pointEstimatePred = varargin{1};
    if nargin > 4
        StartInd =  max(1,varargin{2});
        EndInd   = min(EndInd,varargin{3});
        if StartInd > EndInd
            fprintf('mod_uq_plotSeriesPred_1: StartInd = %g > EndInd = %g \n',...
            StartInd,EndInd);
            error('This relationship is not allowed');
        end % if 
    end % else
end % else 
%% End of block block with modificaation  by O. Klein


% Set custom colors
priorColor = [0.5, 0.7, 1.0];
postColor  = [0.0, 0.2, 0.6];
% O. Klein 2 lines added below
postModColor = [1, 170/255, 1];
postModPEplusColor = [0, 0.5, 1];
%% Create the plot
uq_figure  % Open a new figure

% Use violin plot for histogram of data series
switch plotType
    
    case 'priorPost'
        % Both prior and posterior predictive runs
        %% O. Klein 2 lines modified below
        priorRuns = Sample.PriorPred(:,StartInd:EndInd);
        postRuns = Sample.PostPred(:,StartInd:EndInd);
        % Plot prior predictive runs
        %% O. Klein 1 line modified below
        priorPlot = uq_violinplot(StartInd:EndInd,priorRuns, 'FaceColor', priorColor);
        hold on
        % Plot posterior predictive runs
        %% O. Klein 1 line modified below
        postPlot = uq_violinplot(StartInd:EndInd,postRuns, 'FaceColor', postColor);
        % Update legend target
        legendObj = [priorPlot(1) postPlot(1)];
        legendName = {'prior predictive', 'posterior predictive'};
    case 'prior'
        % Only prior prior predictive runs
        %% O. Klein 1 line modified below
        priorRuns = Sample.PriorPred(:,StartInd:EndInd);
        % Plot prior predictive runs
        %% O. Klein 1 line modified below
        priorPlot = uq_violinplot(StartInd:EndInd,priorRuns, 'FaceColor', priorColor);
        % Update legend target
        legendObj = priorPlot(1);
        legendName = {'prior predictive'};
        hold on
    case 'post'
        % Only posterior predictive runs
        %% O. Klein 1 line modified below
        postRuns = Sample.PostPred(:,StartInd:EndInd);
        % Plot posterior predictive runs
        %% O. Klein 1 line modified below
        postPlot = uq_violinplot(StartInd:EndInd,postRuns, 'FaceColor', postColor);
        % Update legend target
        legendObj = postPlot(1);
        legendName = {'posterior predictive'};
        hold on
        %% O. Klein added new plot version and a consistency check
        %% Start of block with modifications
    case 'postModel'
        % Both posterior predictive and model at samples runs
        postRuns = Sample.PostPred(:,StartInd:EndInd);
        ModelRuns = Sample.Post(:,StartInd:EndInd);
        
 % Plot posterior predictive runs combined with density for model at
 % samples 
        postPlot = uq_violinplot(StartInd:EndInd,postRuns, 'FaceColor', postColor);
        hold on
        % plot model output at sample
        ModelPlot = uq_violinplot(StartInd:EndInd,ModelRuns, 'FaceColor', postModColor);
        hold on
        
        % Update legend target
        legendObj = [postPlot(1) ModelPlot(1)];
        legendName = {'posterior predictive','model at posterior'};
    case 'model'
        % Both posterior predictive and model at samples runs
        ModelRuns = Sample.Post(:,StartInd:EndInd);
        
        
        % plot model output at sample
        ModelPlot = uq_violinplot(StartInd:EndInd,ModelRuns, 'FaceColor', postModColor);
        hold on
        
        % Update legend target
        legendObj = [ModelPlot(1)];
        legendName = {'model at posterior'};
        
        
    case 'Mod(PE)+noise'
        % Both posterior predictive and model at samples runs
        if pointEstimatePredFlag
            % assign model runs to current Data based on MOMap
            pointEstimatePredCurr = zeros(size(pointEstimatePred(1).Out,1),size(Data.y,2));
            % loop over data points
            for ii = 1:size(pointEstimatePredCurr,2)
                currModel = Data.MOMap(1,ii);
                currOut = Data.MOMap(2,ii);
                pointEstimatePredCurr(:,ii) = pointEstimatePred(currModel).Out(:,currOut);
            end
            ModelRuns = Sample.Post(:,StartInd:EndInd);
            postRuns =  Sample.PostPred(:,StartInd:EndInd);
            usedNoise = postRuns -ModelRuns;
            Model_at_PE_plus_noise = pointEstimatePredCurr(1,StartInd:EndInd) + usedNoise;
            
            % plot model output at sample
            Model_at_PE_plusPlot = uq_violinplot(StartInd:EndInd,...
                Model_at_PE_plus_noise, 'FaceColor',postModPEplusColor);
            hold on
            
            % Update legend target
            legendObj = [Model_at_PE_plusPlot(1)];
            legendName = { ['(model at ',lower(pointEstimatePred(1).Type),') +noise']};
        else
            close;
            error('mod_UQLab_plotSeriesPred: Plot option "Mod(PE)+noise" can only be used if point estimate is known');
        end % else to  if    if pointEstimatePredFlag
          case 'postMod(PE)+noise'
        % Plotposterior predictive  density for (model at point estimate ) + noise
        if pointEstimatePredFlag
            % assign model runs to current Data based on MOMap
            pointEstimatePredCurr = zeros(size(pointEstimatePred(1).Out,1),size(Data.y,2));
            % loop over data points
            for ii = 1:size(pointEstimatePredCurr,2)
                currModel = Data.MOMap(1,ii);
                currOut = Data.MOMap(2,ii);
                pointEstimatePredCurr(:,ii) = pointEstimatePred(currModel).Out(:,currOut);
            end
            %% O. Klein 2 lines modified below
            ModelRuns = Sample.Post(:,StartInd:EndInd);
            postRuns =  Sample.PostPred(:,StartInd:EndInd);

            Model_at_PE_plus_noise = pointEstimatePredCurr(1,StartInd:EndInd) + postRuns - ModelRuns;

            % plot model output at sample
                    postPlot = uq_violinplot(StartInd:EndInd,postRuns, 'FaceColor', postColor);

            hold on
                Model_at_PE_plusPlot = uq_violinplot(StartInd:EndInd,...
                Model_at_PE_plus_noise,'FaceColor',postModPEplusColor);
            
            % Update legend target
            legendObj = [postPlot(1) Model_at_PE_plusPlot(1)];
           legendName = {'posterior predictive',...
                ['(model at ',lower(pointEstimatePred(1).Type),') +noise']};
          
        else
        error('mod_UQLab_plotSeriesPred: Plot option "Mod(PE)+noise" can only be used if point estimate is known');
        end % else to  if    if pointEstimatePredFlag
    otherwise
        fprintf('mod_UQLab_plotSeriesPred:\n')
        fprintf('Addmissible values for  plotType are: ''priorPost'', ''prior'',''post'',\n')'
        fprintf(' ''postModel'',  ''Mod(PE)+noise'', and ''postMod(PE)+noise''.');
        error('The used value ``%s`` is not valid,',plotType)
  
    %% End of block with modifications by O.Klein
 
end       




% Use scatter plot for observed data
	%% O. Klein 2 linss  modified below
DataY = Data.y;
xDummy = StartInd:EndInd;
for ii = xDummy
	%% O. Klein 1 line  modified below
    yCurr = DataY(:,ii);
    for jj = 1:numel(yCurr) 
        dataPlot = scatter(ii, yCurr(jj), 100, 'gx'); 
    end
end

% Update legend target
legendObj = [legendObj dataPlot(1)];
legendName{end+1} = 'data';

% Add point estimate
if pointEstimatePredFlag
    % assign model runs to current Data based on MOMap
    pointEstimatePredCurr = zeros(size(pointEstimatePred(1).Out,1),size(Data.y,2));
    % loop over data points
    for ii = 1:size(pointEstimatePredCurr,2)
        currModel = Data.MOMap(1,ii);
        currOut = Data.MOMap(2,ii);
        pointEstimatePredCurr(:,ii) = pointEstimatePred(currModel).Out(:,currOut);
    end
        
	%% O. Klein 1 line  modified below
    pointEstimatePlot = scatter(xDummy, pointEstimatePredCurr(:,xDummy), 100, 'r+');
    % Update legend target
    legendObj = [legendObj,pointEstimatePlot(1)];
    % O. Klein change legend for model outputt at point Estimate
    legendName{end+1} = ['model at ',lower(pointEstimatePred(1).Type)];
end
hold off

% Add legend
uq_legend(legendObj,legendName)

% Add axes labels
ylabel('$\mathcal{Y}$', 'Interpreter', 'LaTeX')
xlabel('$\mathrm{Data\,index\,}$(-)', 'Interpreter', 'LaTeX')

% Set xtick labels (maximum 10)
nTicks = 10;
if length(xDummy) > nTicks
    xDummy = ceil(linspace(1, length(xDummy), nTicks))-1+StartInd;
end
set(gca, 'XTick', xDummy)

% If ticklabel interpreter can be set, update labels to latex
labels = cell(length(xDummy),1);
if isprop(gca,'TickLabelInterpreter')
    for ii = 1:length(xDummy)
        labels{ii} = sprintf('$\\mathrm{y_{%u}}$',xDummy(ii));
    end
else
    for ii = 1:length(xDummy)
        labels{ii} = sprintf('y%d',xDummy(ii));
    end
end
set(gca, 'XTickLabel', labels)

end
