function [Fx, Fy] = pressureForce(p, th, b, r)
%This function evaluates force resulting from pressure field around the wheel
%INPUT
% th: row vector containing angular portions with a n assigned pressure
%     value --> defined fro 0 to 2pi [rad]
% p : row vector containing pressures [Pa] (same length of vector th)
% r : radius of circumference arcs [mm]
% b : gear thickness [mm]
%OUTPUT
% Fx, Fy : force components in XY reference system (driver gear)

b = b*0.001;  %[mm] -> [m]
r = r*0.001;  %[mm] -> [m]
f = 2*p.*sin(th/2)*b*r;  %vector containing force magnitudes resulting from
                         %pressure integration over each circumference arc
phi = zeros(1, length(th)); 
                         %phi is a vector containing force directions 
                         %resulting from pressure integration over each
                         %circumference arc
phi(1) = th(1)/2;
for ii = 2 : length(th)
    phi(ii) = phi(ii-1) + th(ii-1)/2 + th(ii)/2;
end

fx = -f.*cos(phi); fy = -f.*sin(phi); %vectors containing XY force coordinates
Fx = sum(fx); Fy = sum(fy); %resulting pressure Force XY components

end

