function [Prob] = Solver(Prob, SaveType);

%
% Make sure the image directory exists
%
if ~isdir('../images/AngleAssignments')
   !mkdir ../images/AngleAssignments
end

%
% Get Informed PDFs for several solvers
%
% PDF-Types
%    0 - barrier
%    1 - primal simplex
%    2 - dual simplex
%    3 - Weakly Informed Set Cover - scaled by dose rates
%    4 - Weakly Informed Set Cover - unscaled by dose rates
%    5 - Weak Scoring
%
for j = 0:5

%
% We calculate the Balanced PDF from the barrier (interior point) sol'n
%
if j == 0
   Prob.SolverType = j;
   Prob = Solver(Prob, 'cplex');
   Prob.BalancedSoln = BalancedPDF(Prob);
   fprintf('RAD -> Successfully solved for SolverType %d \n', j)
   Prob.InteriorAvg = zeros(1,Prob.NumAngles);
   Prob.BalancedAvg = zeros(1,Prob.NumAngles);
   Prob.InteriorMax = zeros(1,Prob.NumAngles);
   Prob.BalancedMax = zeros(1,Prob.NumAngles);
   Prob.EntropyAvgBalanced = zeros(1, Prob.NumAngles);
   ColInd = 1;
   for a = 1:Prob.NumAngles
      Prob.InteriorAvg(a) = ...
         sum(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      Prob.BalancedAvg(a) = ...
         sum(Prob.BalancedSoln( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      Prob.InteriorMax(a) = ...
         max(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      Prob.BalancedMax(a) = ...
         max(Prob.BalancedSoln( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      ColInd = ColInd + Prob.NumPencils(a);
   end
   ColInd = 1;
   for a = 1:Prob.NumAngles
      % Calculate angles Entropy
      Prob.EntropyAvgBalanced(a) = sum(...
         Prob.BalancedSoln(ColInd : ColInd + Prob.NumPencils(a) -1) .* ...
         log(max(Prob.BalancedSoln(ColInd : ColInd + Prob.NumPencils(a) -1), ...
         10^(-10)*ones(Prob.NumPencils(a), 1)) ) );
      ColInd = ColInd + Prob.NumPencils(a);
   end
   % Calculate Entropy PDFs
   Prob.EntropyAvgBalanced = (1/max(Prob.EntropyAvgBalanced)) * ...
      (Prob.EntropyAvgBalanced - min(Prob.EntropyAvgBalanced));
   %
   % Plot and save these PDfs
   %
   figure
   plot(1:Prob.NumAngles, Prob.InteriorAvg);
   saveas(gcf,['../images/AngleAssignments/InteriorAvg'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.InteriorMax);
   saveas(gcf,['../images/AngleAssignments/InteriorMax'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.BalancedAvg);
   saveas(gcf,['../images/AngleAssignments/BalancedAvg'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.BalancedMax);
   saveas(gcf,['../images/AngleAssignments/BalancedMax'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.EntropyAvgBalanced);
   saveas(gcf,['../images/AngleAssignments/EntropyAvgBalanced'], SaveType)

%
% Find the PDFs for the primal simplex
%
elseif j == 1
   Prob.SolverType = j;
   Prob = Solver(Prob, 'cplex');
   fprintf('RAD -> Successfully solved for SolverType %d \n', j)
   Prob.PrimalAvg = zeros(1,Prob.NumAngles);
   Prob.PrimalMax = zeros(1,Prob.NumAngles);
   ColInd = 1;
   for a = 1:Prob.NumAngles
      Prob.PrimalAvg(a) = ...
         sum(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1));
      Prob.PrimalMax(a) = ...
         max(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1));
      ColInd = ColInd + Prob.NumPencils(a);
   end
   figure
   plot(1:Prob.NumAngles, Prob.PrimalAvg);
   saveas(gcf,['../images/AngleAssignments/PrimalAvg'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.PrimalMax);
   saveas(gcf,['../images/AngleAssignments/PrimalMax'], SaveType)
%
% Find the PDFs for the Dual simplex
%
elseif j == 2
   Prob.SolverType = j;
   Prob = Solver(Prob, 'cplex');
   fprintf('RAD -> Successfully solved for SolverType %d \n', j)
   Prob.DualAvg = zeros(1,Prob.NumAngles);
   Prob.DualMax = zeros(1,Prob.NumAngles);
   ColInd = 1;
   for a = 1:Prob.NumAngles
      Prob.DualAvg(a) = ...
         sum(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      Prob.DualMax(a) = ...
         max(Prob.ThePlan( ColInd : ColInd + Prob.NumPencils(a) - 1 ));
      ColInd = ColInd + Prob.NumPencils(a);
   end
   figure
   plot(1:Prob.NumAngles, Prob.DualAvg);
   saveas(gcf,['../images/AngleAssignments/DualAvg'], SaveType)
   figure
   plot(1:Prob.NumAngles, Prob.DualMax);
   saveas(gcf,['../images/AngleAssignments/DualMax'], SaveType)
%
% Find the PDF for the Scaled weakly informed set cover technique.
%
elseif j == 3
   Prob.SetCover1 = SetCovPDF(Prob,'ScaledWeak',[]);
   figure
   plot(1:Prob.NumAngles, Prob.SetCover1);
   saveas(gcf,['../images/AngleAssignments/ScaledWeakSetCover'], SaveType)
   fprintf('RAD -> Successfully Calculated Scaled Weak SetCover PDF \n')
%
% Find the PDF for the Unscaled weakly informed set cover technique.
%
elseif j == 4
   Prob.SetCover2 = SetCovPDF(Prob, 'UnScaledWeak', []);
   figure
   plot(1:Prob.NumAngles, Prob.SetCover2);
   saveas(gcf,['../images/AngleAssignments/UnscaledWeakSetCover'], SaveType)
   fprintf('RAD -> Successfully Calculated Unscaled Weak SetCover PDF \n')
%
% Find the PDFs for the scoring technique
%
elseif j == 5
   Prob.Scoring = SetCovPDF(Prob, 'Scoring', []);
   figure
   plot(1:Prob.NumAngles, Prob.Scoring);
   saveas(gcf,['../images/AngleAssignments/Scoring'], SaveType)
   fprintf('RAD -> Successfully Calculated Scoring PDF \n')
else
   fprintf('RAD -> WARNING in PDFmaker: Unknown Solver Type \n')
end

%
% End for loop for different solvers
%   
end
