function poweran(filename, pf_angle)
	% Power analyser
	%
	% Analyse and plot the instantaneous, real, and reactive
	% power of a waveform given by a tab-separated value
	% (TSV) file.
	%
	% Required measurements to input:
	% power factor angle (degree)
	%
	% The following column headers for the TSV are assumed:
	% time, v1, v2, v3, i1, i2, i3
	%
	% Time measurements are assumed to be in milliseconds.
	% All voltages and currents are assumed to be given
	% in Volts and Amperes, respectively.

	data = readtable(filename, FileType="text", Delimiter="\t");

	L = length(data.time);

	t = linspace(0, data.time(end), L);
	p = data.v1(1:L) .* data.i1(1:L);

	% Sum over a full period to get the RMS values.
	n = length(data.time(data.time < 1e3/60));
	vrms = sqrt(sum(data.v1(data.time < 1e3/60).^2)/n);
	irms = sqrt(sum(data.i1(data.time < 1e3/60).^2)/n);

	% These are constant in time so we extrapolate with two points
	% to save compute time.
	real_power = ones(2) * vrms * irms * cos(deg2rad(pf_angle));
	reactive_power = ones(2) * vrms * irms * sin(deg2rad(pf_angle));

	display(vrms);
	display(irms);

	figure(1);

	subplot(3, 1, 1);
	plot(t, p, "r", LineWidth=2);
	title("Instantaneous Power versus Time");
	xlabel("Time, $t$ (ms)", Interpreter="latex");
	ylabel("Power, $p(t)$ (W)", Interpreter="latex");

	subplot(3, 1, 2);
	plot([0 data.time(end)], real_power, "g", LineWidth=2);
	title("Real Power versus Time");
	xlabel("Time, $t$ (ms)", Interpreter="latex");
	ylabel("Real Power, $P$ (W)", Interpreter="latex");

	subplot(3, 1, 3);
	plot([0 data.time(end)], reactive_power, "b", LineWidth=2);
	title("Reactive Power versus Time");
	xlabel("Time, $t$ (ms)", Interpreter="latex");
	ylabel("Reactive Power, $Q$ (VAR)", Interpreter="latex");
end