function it_pounds, it_pct, nic_pounds, nic_pct, total_pounds, total_pct, takehome_pounds = tax(gross, tax_rates, tax_widths, nic_rates, nic_widths, plot_all)
%TAX: Compute the pounds and percentage cost of taxable pay of income tax
%and National Insurance. (Or just tax in places without NI).
%
%[it_pounds it_pct nic_pounds nic_pct total_pounds total_pct takehome_pounds] = tax(gross, tax_rates, tax_widths, nic_rates, nic_widths, plot_all)
%
%Inputs:
%   gross:      Vector of taxable pay amounts to compute all outputs for.
%   tax_rates:  Percentages of tax rates matching entries in tax_widths.
%   tax_widths: Width in pounds of each tax band, ended by Inf. NOT the
%   income levels at which they begin to apply. For example, for 2007-08:
%
%   Personal allowance, or 0% band          = £5225.
%   10% band = First £2230 taxable          = £2230 in width.
%   22% band = £2230 - £34600 taxable       = £32370 in width.
%   40% band = £34600 - (unlimited) taxable = £Inf in width.
%
%   Therefore, tax_rates and tax_widths are set up as follows:
%
%   tax_rates  = [ 0    10    22   40];
%   tax_widths = [5225 2230 32370 Inf];
%
%   nic_rates:  As tax_rates. (Omit if no NIC-equivalent needed).
%   nic_widths: As tax_widths. (Omit if no NIC-equivalent needed).
%   plot_all :  Single boolean true/false whether to draw the output plot.
%
%Outputs:
%   These all have obvious names, and each produce a vector with entries
%   corresponding to the entries of input 'gross'.
 
for g = 1:length(gross)
    it_pounds(g) = tax_calc(gross(g), tax_rates, tax_widths);
    it_pct(g) = 100*it_pounds(g)/gross(g);
 
    nic_pounds(g) = tax_calc(gross(g), nic_rates, nic_widths);
    nic_pct(g) = 100*nic_pounds(g)/gross(g);
 
    total_pounds(g) = it_pounds(g) + nic_pounds(g);
    total_pct(g) = 100*total_pounds(g)/gross(g);
 
    takehome_pounds(g) = gross(g) - total_pounds(g);
end
 
if exist('plot_all', 'var') && plot_all
 
    %Plots for amounts in pounds
    figure;
    plot(gross, it_pounds, 'b', gross, nic_pounds, 'r', gross, total_pounds, 'k', gross, takehome_pounds, 'g');
    set(get(gca, 'Children'), 'LineWidth',2);
    set(gca, 'XMinorTick', 'on', 'YMinorTick', 'on');
    legend('Income tax', 'National insurance', 'Combined', 'Take home');
    xlabel('Annual taxable pay, GBP');
    ylabel('Annual amount, GBP');
    grid on;
 
    %Plots for percentages
    figure;
    plot(gross, it_pct, 'b', gross, nic_pct, 'r', gross, total_pct, 'k');
    set(get(gca, 'Children'), 'LineWidth',2);
    set(gca, 'XMinorTick', 'on', 'YMinorTick', 'on');
    legend('Income tax', 'National insurance', 'Combined');
    xlabel('Annual taxable pay, GBP');
    ylabel('Percentage of taxable pay');
    grid on;
end
 
%%%%Calculation engine%%%%
function total_tax = tax_calc(gross, rates, widths)
 
tax_band = 1;
income_left = gross;
total_tax = 0;
 
while tax_band <= length(widths) && income_left > 0
    taxable_this_band = min(widths(tax_band) income_left);
    total_tax = total_tax + taxable_this_band*rates(tax_band)/100;
    income_left = income_left - taxable_this_band;
    tax_band = tax_band + 1;
end

tax_rates  = 0 20 40;
tax_widths = 6035 34800 Inf; %until changes announced [[13 May]] [[2008]], this was [5435 36000 Inf].
nic_rates  = 0 11 1;
nic_widths = 105*52 (770-105)*52 Inf;