# Puka - Software for Detection of Breaths in Strain Gauge Recordings 1.0.0

File: <base>/puka/matlabScripts/newPT.m (2,807 bytes)
``````function [P,T,th,Qd] = newPT(Qraw, factor, onsetTime, endTime);
% have the respiration signal in a column vector called Q
% call by >> [P,T,th,Qd] = newPT(Q);
% code adapted from Todd & Andrews. The Identification of Peaks in Physiological
% Signals. Computers and Biomedical Research 32, 322-335 (1999).

Qd = decimate(Qraw, 5);  % downsample the signal

d = 1;      % variable to show if signal unknown (1), trough-to-peak (2),
% or peak-to-trough (3)
a = 1;      % index of maximal element since last trough
b = 1;      % index of minimal element since last peak
S = [];     % indices of maximal elements since last trough if signal rising (d = 2)
% or minimal elements since last peak if signal falling (d = 3)
P = [];     % peak elements
T = [];     % trough elements
% factor is the number to use when making the threshold; default is 0.1

th = abs(prctile(Qd,75) - prctile(Qd,25)) * factor;
% th = abs(max(Qd) - min(Qd)) * .5;  % threshold

[n] = max(size(Qd));  % get size of the array so can do the for loop
for i = 1:n
if d == 1
if Qd(a) >= Qd(i) + th
d = 3;
elseif Qd(i) >= Qd(b) + th
d = 2;
end;
if Qd(a) < Qd(i)
a = i;
elseif Qd(i) < Qd(b)
b = i;
end;
S = i;
elseif d == 2  % signal rising, trough-to-peak
if Qd(a) < Qd(i) % still rising
S = i;
a = i;
elseif Qd(a) == Qd(i)
S = [S,i];
elseif Qd(a) >= Qd(i) + th
P = [P,S]; S = i;
b = i; d = 3;
end;
elseif d == 3  % signal falling, peak-to-trough
if Qd(i) <= Qd(b)
S = i;
b = i;
elseif Qd(b) == Qd(i)
S = [S,i];
elseif Qd(i) >= Qd(b) + th
T = [T,S]; S = i;
a = i; d = 2;
end;
end;
end;

% go through P and T, removing marks for peaks/troughs before the start point
% and after the end point
goodP = [];  goodT = [];
onsetTime = onsetTime/5; endTime = endTime/5;  % to match decimation
n = length(P);
for i = 1:n
if P(i) > onsetTime
if P(i) < endTime
goodP = [goodP, P(i)];   % only add in peaks bigger than onsetTime
end;
end;
end;
n = length(T);
for i = 1:n
if T(i) > onsetTime
if T(i) < endTime
goodT = [goodT, T(i)];   % only add in troughs bigger than onsetTime
end;
end;
end;

P = []; T = []; P = goodP; T = goodT;  % reset arrays for return values

plot(Qd, 'm'); whitebg([.9 .9 .9]);   % set background color to gray
hold on;  % plot the signal and peaks
plot(P, Qd(P), 'ob', 'MarkerSize',10);
plot(T, Qd(T), 'ob', 'MarkerSize',10);
hold off;``````