Plots the schedule's plans, cross sections, production, grades
--------------------------------------------------------------------------
Date: Aug 06, 2008, ver01
By: Hooman Askari
--------------------------------------------------------------------------
Inputs
Blocks:
a struct containing Blocks(numBlocks) elements with the field
correspoiding to the input file format.
Result:
a struct holding the CPLEX solution. the solution is in the
field x_k
numBlocks:
a 1*1 vector containing the total number of blocks
numOfPeriods:
a 1*1 vector containing the number of periods
normC:
the norm of the objective function coefficients
xdim; ydim; zdim:
small block size in x, y, and z direction
--------------------------------------------------------------------------
--------------------------------------------------------------------------
outputs
Plots (per period) the:
Total tonnage mined ; Ore tonnage;
Waste tonnage; Tonnage of Ore and Waste;
Cash Flow; Average grade of ore and contaminants;
Period: Period(1:numBlocks, 1)= holds the period that block i has
been last extracted.
scheduledBlocks: save the results as a struct with the format of Blocks.
with adding the array Blocks.periods(1:2,:) to the
Blocks struct. Blocks.periods(1 ,:) = periods that
each block is extracted in Blocks.periods(2 ,:) =
portion of the block that is extracted in the
corrosponding period.
graphs: is a struct that saves the the attributes of the plotted
variables per year
--------------------------------------------------------------------------
--------------------------------------------------------------------------
Pseudo Code
--------------------------------------------------------------------------0001 % Plots the schedule's plans, cross sections, production, grades 0002 %-------------------------------------------------------------------------- 0003 % Date: Aug 06, 2008, ver01 0004 % By: Hooman Askari 0005 %-------------------------------------------------------------------------- 0006 %Inputs 0007 % Blocks: 0008 % a struct containing Blocks(numBlocks) elements with the field 0009 % correspoiding to the input file format. 0010 % Result: 0011 % a struct holding the CPLEX solution. the solution is in the 0012 % field x_k 0013 % numBlocks: 0014 % a 1*1 vector containing the total number of blocks 0015 % numOfPeriods: 0016 % a 1*1 vector containing the number of periods 0017 % normC: 0018 % the norm of the objective function coefficients 0019 % xdim; ydim; zdim: 0020 % small block size in x, y, and z direction 0021 %-------------------------------------------------------------------------- 0022 %-------------------------------------------------------------------------- 0023 % outputs 0024 % 0025 % Plots (per period) the: 0026 % Total tonnage mined ; Ore tonnage; 0027 % Waste tonnage; Tonnage of Ore and Waste; 0028 % Cash Flow; Average grade of ore and contaminants; 0029 % 0030 % Period: Period(1:numBlocks, 1)= holds the period that block i has 0031 % been last extracted. 0032 % 0033 % scheduledBlocks: save the results as a struct with the format of Blocks. 0034 % with adding the array Blocks.periods(1:2,:) to the 0035 % Blocks struct. Blocks.periods(1 ,:) = periods that 0036 % each block is extracted in Blocks.periods(2 ,:) = 0037 % portion of the block that is extracted in the 0038 % corrosponding period. 0039 % graphs: is a struct that saves the the attributes of the plotted 0040 % variables per year 0041 %-------------------------------------------------------------------------- 0042 % 0043 %-------------------------------------------------------------------------- 0044 % Pseudo Code 0045 % 0046 %-------------------------------------------------------------------------- 0047 function output_plot_cut(Blocks, Result, numOfPeriods, numBlocks,... 0048 numCuts, normC, xdim, ydim, zdim) 0049 0050 % steps that will take to plot 0051 step = 3; 0052 0053 % xdim = dimension of block in the X direction, Easting; 0054 % ydim = dimension of block in the Y direction, Northing; 0055 % zdim = dimension of block in the Z direction, Elevation; 0056 0057 % algorithm used by sovler 0058 solverAlgorithm = Result.SolverAlgorithm 0059 0060 % CPU time used in seconds 0061 cpuTime = Result.CPUtime 0062 0063 % number of major iterations 0064 % numOfIterations = Result.iter 0065 0066 % Number of constraint evaluations needed. 0067 constraintEvalualtions = Result.ConstrEv 0068 0069 % discounted cash flow at the optimum schedule 0070 display('Discounted Cash Flow in Million Dollars') 0071 DCF = cur2str((Result.f_k*-1* normC)/1000000,2) 0072 0073 %-------------------------------------------------------------------------- 0074 % binSchedule [numBlocks, numOfPeriods]= binary integer matrix holding the 0075 % schedule 0076 % each column contains the corresponding binary schedule of that period. 0077 % e.g. second column holds the second period schedule 0078 0079 % b(K*T) = Result.x_k( 1: K*T )= holds the binary integer variable 0080 % y(K*T) = Result.x_k( [end-(K*T)+1] : end )= holds the continuous cut 0081 % variable 0082 % x(N*T) = Result.x_k( [(K*T)+1] : [end-(K*T)] )= holds the continuous block 0083 % variable 0084 0085 x = Result.x_k( (numCuts*numOfPeriods+1) : (end-(numCuts*numOfPeriods))); 0086 0087 binSchedule = reshape(x, numBlocks, numOfPeriods); 0088 0089 % for iPeriods = 1: numOfPeriods 0090 % binSchedule(:,iPeriods) = Result.x_k((numCuts * (iPeriods))+ 1:... 0091 % (numCuts * (iPeriods-1))+ numBlocks ,1); 0092 % end 0093 0094 [m n]=size(binSchedule); 0095 binWasteSchedule = zeros(m, n); 0096 binOreSchedule = zeros(m, n); 0097 0098 % constructing matricies with the same format of binSchedule used for plots 0099 for iPeriods = 1: numOfPeriods 0100 for iBlocks = 1: numBlocks 0101 blockTonnage(iBlocks,iPeriods) = Blocks(iBlocks).BlockTonnage; 0102 EBV(iBlocks,iPeriods) = Blocks(iBlocks).EBV; 0103 0104 if (EBV(iBlocks,iPeriods)<=0 && binSchedule(iBlocks,iPeriods)~= 0) 0105 % binWasteSchedule(iBlocks,iPeriods) = 1; 0106 binWasteSchedule(iBlocks,iPeriods) = binSchedule(iBlocks,iPeriods); 0107 elseif (EBV(iBlocks,iPeriods)> 0 && binSchedule(iBlocks,iPeriods)~= 0) 0108 %binOreSchedule(iBlocks,iPeriods) = 1; 0109 binOreSchedule(iBlocks,iPeriods) = binSchedule(iBlocks,iPeriods); 0110 end 0111 0112 MWT(iBlocks,iPeriods) = Blocks(iBlocks).MWT; 0113 S(iBlocks,iPeriods) = Blocks(iBlocks).S; 0114 P(iBlocks,iPeriods) = Blocks(iBlocks).P; 0115 gradeMWT(iBlocks,iPeriods) = Blocks(iBlocks).gradeMWT; 0116 gradeS(iBlocks,iPeriods) = Blocks(iBlocks).gradeS; 0117 gradeP(iBlocks,iPeriods) = Blocks(iBlocks).gradeP; 0118 OreTonnes(iBlocks,iPeriods) = Blocks(iBlocks).OreTonnes; 0119 WasteTonnes(iBlocks,iPeriods) = Blocks(iBlocks).WasteTonnes; 0120 OreValue(iBlocks,iPeriods) = Blocks(iBlocks).oreValue; 0121 WasteCosts(iBlocks,iPeriods) = Blocks(iBlocks).miningCost; 0122 0123 end 0124 end 0125 0126 tonnageMined = binSchedule .* blockTonnage; 0127 cashFlow = binSchedule .* EBV; 0128 % ore send to mill 0129 % productTonnage = binSchedule .* MWT; 0130 oreTonnage = binSchedule .* OreTonnes; 0131 wasteTonnage = binSchedule .* WasteTonnes; 0132 sTonnage = binSchedule .* S; 0133 pTonnage = binSchedule .* P; 0134 0135 % oreGrade = binSchedule .*gradeMWT; 0136 % sGrade = binSchedule .* gradeS; 0137 % pGrade = binSchedule .* gradeP; 0138 0139 oreGrade = gradeMWT; 0140 sGrade = gradeS; 0141 pGrade = gradeP; 0142 0143 % Plotting grades, etc----------------------------------------------------- 0144 % createBar(data, xLabel, yLabel, Title) 0145 0146 createBar(sum(tonnageMined)/1000000,... 0147 'Period',... 0148 'Tonnage Mined(MTonne)',... 0149 'Total Tonnage Mined'); 0150 0151 createBar(sum(oreTonnage)/1000000,... 0152 'Period',... 0153 'Ore Tonnage (MTonne)',... 0154 'Ore Tonnage'); 0155 0156 createBar(sum(wasteTonnage)/1000000,... 0157 'Period',... 0158 'Waste Tonnage (MTonne)',... 0159 'Waste Tonnage'); 0160 0161 Y = [sum(oreTonnage)/1000000; sum(wasteTonnage)/1000000]'; 0162 createBar(Y,... 0163 'Period',... 0164 'Ore and waste tonnage (MTonne)',... 0165 'Tonnage of Ore and Waste'); 0166 0167 0168 createBar(sum(cashFlow)/1000000,... 0169 'Period',... 0170 'cash flow ($M)',... 0171 'Cash Flow'); 0172 % %---------------------------------------------------------------------- 0173 % UNCOMMENT this section to plot high resoultion Figures 0174 % Plotting the same figures as above with high resolution as TIF files 0175 % 0176 % createBar(sum(tonnageMined)/1000000,... 0177 % 'Period',... 0178 % 'Tonnage Mined(MTonne)',... 0179 % 'Total Tonnage Mined'); 0180 % 0181 % print -dtiff -r150 im.tiff 0182 % im = imread( 'im.tiff', 'tiff' ); 0183 % [im_hatch,colorlist] = applyhatch_pluscolor(im,'\-x.',1,[],[],100,3,3); 0184 % imwrite(im_hatch,'im_hatch.tiff','tiff') 0185 % 0186 % createBar(sum(oreTonnage)/1000000,... 0187 % 'Period',... 0188 % 'Ore Tonnage (MTonne)',... 0189 % 'Ore Tonnage'); 0190 % 0191 % print -dtiff -r150 im.tiff 0192 % im = imread( 'im.tiff', 'tiff' ); 0193 % [im_hatch,colorlist] = applyhatch_pluscolor(im,'\-x.',1,[],[],100,3,3); 0194 % imwrite(im_hatch,'im_hatch.tiff','tiff') 0195 % 0196 % createBar(sum(wasteTonnage)/1000000,... 0197 % 'Period',... 0198 % 'Waste Tonnage (MTonne)',... 0199 % 'Waste Tonnage'); 0200 % 0201 % print -dtiff -r150 im.tiff 0202 % im = imread( 'im.tiff', 'tiff' ); 0203 % [im_hatch,colorlist] = applyhatch_pluscolor(im,'\-x.',1,[],[],100,3,3); 0204 % imwrite(im_hatch,'im_hatch.tiff','tiff') 0205 % 0206 % Y = [sum(oreTonnage)/1000000; sum(wasteTonnage)/1000000]'; 0207 % createBar(Y,... 0208 % 'Period',... 0209 % 'Ore and waste tonnage (MTonne)',... 0210 % 'Tonnage of Ore and Waste'); 0211 % 0212 % print -dtiff -r150 im.tiff 0213 % im = imread( 'im.tiff', 'tiff' ); 0214 % [im_hatch,colorlist] = applyhatch_pluscolor(im,'\-x.',1,[],[],100,3,3); 0215 % imwrite(im_hatch,'im_hatch.tiff','tiff') 0216 % 0217 % createBar(sum(cashFlow)/1000000,... 0218 % 'Period',... 0219 % 'cash flow (M$)',... 0220 % 'Cash Flow'); 0221 % 0222 % print -dtiff -r150 im.tiff 0223 % im = imread( 'im.tiff', 'tiff' ); 0224 % [im_hatch,colorlist] = applyhatch_pluscolor(im,'\-x.',1,[],[],100,3,3); 0225 % imwrite(im_hatch,'im_hatch.tiff','tiff') 0226 % %------------------------------------------------------------------------ 0227 0228 for iPeriods = 1: numOfPeriods 0229 0230 % SumProduct of ore grade and tonnage to calculate average grade 0231 sumProductGradeOreTonnage(iPeriods) = ... 0232 sum(oreGrade(:,iPeriods).* oreTonnage(:,iPeriods)); 0233 % sum of ore tonnage 0234 sumOreTonnage(iPeriods) = sum(oreTonnage(:,iPeriods)); 0235 0236 sumProductSOreTonnage(iPeriods) = ... 0237 sum(sGrade(:,iPeriods).* sTonnage(:,iPeriods)); 0238 sumSTonnage(iPeriods) = sum(sTonnage(:,iPeriods)); 0239 0240 sumProductPOreTonnage(iPeriods) = ... 0241 sum(pGrade(:,iPeriods).* pTonnage(:,iPeriods)); 0242 sumPTonnage(iPeriods) = sum(pTonnage(:,iPeriods)); 0243 0244 end 0245 %-------------------------------------------------------------------------- 0246 % createPlot(xData, yData, xLabel, yLabel, Title) 0247 0248 createPlot( 1:numOfPeriods,... 0249 (sumProductGradeOreTonnage./sumOreTonnage)*100,... 0250 'Period',... 0251 ' MWT% ',... 0252 'Average grade MWT% per period'); 0253 0254 createPlot( 1:numOfPeriods,... 0255 (sumProductSOreTonnage./sumSTonnage)*100,... 0256 'Period',... 0257 'S%',... 0258 'Average grade S% per period'); 0259 0260 createPlot( 1:numOfPeriods,... 0261 (sumProductPOreTonnage./sumPTonnage)*100,... 0262 'Period',... 0263 ' P% ',... 0264 'Average grade P per period'); 0265 0266 % struct Grpahs holds the attributes per year------------------------------ 0267 Graphs = struct(... 0268 'oreTonnagePerYr',sum(oreTonnage)/1000000,... 0269 'wasteTonnagePerYr' , sum(wasteTonnage)/1000000,... 0270 'cashFlowPerYr', sum(cashFlow)/1000000,... 0271 'oreGradePerYr', (sumProductGradeOreTonnage./sumOreTonnage)*100,... 0272 'sGradePerYr', (sumProductSOreTonnage./sumSTonnage)*100,... 0273 'pGradePerYr', (sumProductPOreTonnage./sumPTonnage)*100,... 0274 'tonnageMinedPerYr',(sum(tonnageMined)/1000000)); 0275 0276 save graphs Graphs; 0277 0278 for iBlocks = 1: numBlocks 0279 0280 % adds the array Blocks.periods(1:2,:) to the Blocks struct 0281 % Blocks.periods(1 ,:) = periods that each block is extracted in 0282 % Blocks.periods(2 ,:) = portion of the block that is extracted in the 0283 % corrosponding period. 0284 0285 % returns the periods 0286 Blocks(iBlocks).periods(1, :) = find (binSchedule(iBlocks, :)); 0287 0288 % returns the portion of the block extracted in each period 0289 Blocks(iBlocks).periods(2, :) = binSchedule(iBlocks, ... 0290 Blocks(iBlocks).periods(1, :)); 0291 end % end for iBlocks 0292 0293 % save the Blocks with their respective extraction periods. 0294 save scheduledBlocks Blocks 0295 0296 %-------------------------------------------------------------------------- 0297 % finding the maximum and minimum block indices to plot the schedule 0298 x = zeros(numBlocks, 1); y = zeros(numBlocks, 1); z = zeros(numBlocks, 1); 0299 0300 for iBlocks = 1: numBlocks 0301 x(iBlocks) = Blocks(iBlocks).XI; 0302 y(iBlocks) = Blocks(iBlocks).YI; 0303 z(iBlocks) = Blocks(iBlocks).ZI; 0304 end % end for iBlocks 0305 0306 % max_XI holds the maximum index of blocks in X direction 0307 % max_XI_index holds the maximimum general index (iBlocks = 1:numBlocks) 0308 % of block max_XI 0309 [max_XI, max_XI_index] = max(x); [min_XI, min_XI_index] = min(x); 0310 [max_YI, max_YI_index] = max(y); [min_YI ,min_YI_index] = min(y); 0311 [max_ZI, max_ZI_index] = max(z); [min_ZI, min_ZI_index] = min(z); 0312 0313 % coordinates of the origin of the blocks that will be plotted 0314 xCoorMin = Blocks(min_XI_index).X; 0315 yCoorMin = Blocks(min_YI_index).Y; 0316 zCoorMin = Blocks(min_ZI_index).Z; 0317 0318 xCoorMax = Blocks(max_XI_index).X; 0319 yCoorMax = Blocks(max_YI_index).Y; 0320 zCoorMax = Blocks(max_ZI_index).Z; 0321 %-------------------------------------------------------------------------- 0322 % switching the Matlab matrix coordinates to match the whittle coordinates 0323 % check documentation for the reason 0324 max_matrix_XI = max_YI; 0325 max_matrix_YI = max_XI; 0326 max_matrix_ZI = max_ZI; 0327 0328 min_matrix_XI = min_YI; 0329 min_matrix_YI = min_XI; 0330 min_matrix_ZI = min_ZI; 0331 0332 blockModel = (numOfPeriods + 3) * ones( max_matrix_XI,... 0333 max_matrix_YI,... 0334 max_matrix_ZI ); 0335 for iBlocks = 1: numBlocks 0336 0337 % see documentation for this section 0338 % XI location based on Whittle coordinate 0339 % matrix_XI coordinates 0340 0341 % these are for checking delete later on 0342 % Blocks(iBlocks).XI; 0343 % Blocks(iBlocks).YI; 0344 % Blocks(iBlocks).ZI; 0345 0346 matrix_XI = max_matrix_XI - Blocks(iBlocks).YI + 1; 0347 matrix_YI = Blocks(iBlocks).XI; 0348 matrix_ZI = Blocks(iBlocks).ZI; 0349 0350 % this is for checking delete later on 0351 % period = Blocks(iBlocks).periods(1,1) 0352 0353 % For the purpose of visualizing the blocks a decision must be made 0354 % for the blocks that are extracted in different periods. here we are 0355 % consdidering the period that has the largest fraction for 0356 % visualizing this can cause some problems in the order of the blocks 0357 % that are shown! be aware and make sure it is working fine. 0358 0359 [maxVal index] = max(Blocks(iBlocks).periods(2,:)); 0360 0361 blockModel(matrix_XI, matrix_YI, matrix_ZI) = ... 0362 Blocks(iBlocks).periods(1,index); 0363 0364 % % OR the first period that a portion of a block is extracted 0365 % blockModel(matrix_XI, matrix_YI, matrix_ZI) = ... 0366 % Blocks(iBlocks).periods(1,end); 0367 0368 % save the extraction period in Period to export to the excel file 0369 Period(iBlocks) = Blocks(iBlocks).periods(1,1); 0370 0371 end % end for iBlocks 0372 0373 % transpose 0374 % Period(1:numBlocks, 1)= holds the period that block i has been last 0375 % extracted 0376 Period = Period'; 0377 save Period Period; 0378 0379 blockModel = blockModel( (max_YI - max_matrix_XI + 1):... 0380 (max_matrix_XI - min_YI + 1),... 0381 min_matrix_YI: max_matrix_YI,... 0382 min_matrix_ZI: max_matrix_ZI ); 0383 0384 % m in the Y direction - Northing 0385 % n in the X direction - Easting 0386 % r in the Z direction - Eleveation 0387 [m n r] = size(blockModel); 0388 0389 % this is how the block model is setup 0390 % blockModel(Y, X, Z) 0391 0392 % Elevation - Planview 0393 for k = 1:r 0394 temp = blockModel(:,:,k); 0395 temp = squeeze(temp); 0396 figure 0397 % axis([xCoorMin xCoorMax yCoorMin yCoorMax]); 0398 viewmatrix(temp) 0399 xlabel( {'Easting (X-index)'; ['each cell ' num2str(xdim) 'm']},... 0400 'FontName','Times New Roman',... 0401 'FontSize',14); 0402 ylabel( {'Northing (Y-index)'; ['each cell ' num2str(ydim) 'm']},... 0403 'FontName','Times New Roman',... 0404 'FontSize',14); 0405 title(['Plan View - ' num2str(zCoorMin + (k-1)*zdim) 'm'],... 0406 'FontName','Times New Roman',... 0407 'FontSize',14); 0408 colormap(hot) 0409 set(gca,'DataAspectRatio',[ydim xdim 1]) 0410 end 0411 0412 % cross sections looking East on the X axis 0413 for j = 1:step:n 0414 temp = blockModel(:,j,:); 0415 temp = squeeze(temp); 0416 temp = rot90(temp); 0417 figure 0418 viewmatrix(temp) 0419 xlabel({'Northing (Y-index)'; ['each cell ' num2str(ydim) 'm']},... 0420 'FontName','Times New Roman',... 0421 'FontSize',14); 0422 ylabel({'Elevation (Z-index)'; ['each cell ' num2str(zdim) 'm']},... 0423 'FontName','Times New Roman',... 0424 'FontSize',14); 0425 title(['Cross Section Looking East - '... 0426 num2str(xCoorMin + (j-1)*xdim) 'm'],... 0427 'FontName','Times New Roman',... 0428 'FontSize',14); 0429 colormap(hot) 0430 set(gca,'DataAspectRatio',[zdim ydim 1]) 0431 end 0432 0433 % cross sections looking North on the Y axis 0434 for i = 1:step:m 0435 temp = blockModel(i,:,:); 0436 temp = squeeze(temp); 0437 temp = rot90(temp); 0438 figure 0439 viewmatrix(temp) 0440 xlabel({'Easting (X-index)'; ['each cell ' num2str(xdim) 'm']},... 0441 'FontName','Times New Roman',... 0442 'FontSize',14); 0443 ylabel({'Elevation (Z-index)'; ['each cell ' num2str(zdim) 'm']},... 0444 'FontName','Times New Roman',... 0445 'FontSize',14); 0446 title(['Cross Section Looking North - '... 0447 num2str(yCoorMin + (i-1)*ydim) 'm'],... 0448 'FontName','Times New Roman',... 0449 'FontSize',14); 0450 colormap(hot) 0451 set(gca,'DataAspectRatio',[zdim xdim 1]) 0452 end 0453 0454 end