function [angle] = meshing_gears(tcx, tcy, O1, O2, PA, Z)
%% Questa funzione calcola la rotazione della ruota condotta 
%(in senso antiorario) per garantire il contatto all'ingranamento
%Sistema di riferimento centrale O(0,0)
% INPUT: 
% 1. tcx, tcy : profilo del dente (lato sx e dx)
% 2. O1, O2 : posizioni di assi rotazione ruota 1,2 (coordinata x)
% 3. PA : angolo di pressione
% 4. Z : numero di denti 
%OUTPUT: 
% 1. angle, angolo di cui ruotare la ruota condotta
% 2. PLOT di tutte le posizioni angolari fino a matching condotta-motrice

[X1, Y1] = rotateZ(tcx(length(tcx)/2 + 1 : end), ...
                   tcy(length(tcy)/2 + 1 : end), -pi/2);
X1 = X1 + O1;
[X2, Y2] = rotateZ(tcx(length(tcx)/2 + 1 : end), ...
                   tcy(length(tcy)/2 + 1 : end), pi/2);
Xm = 0; Ym = 0;

tol = abs(1e-3 * O1);
%Valutazione del punto di contatto, giacente su linea dei contatti definita
%da angolo di pressione 
for i = 1:length(X1)
    if ( abs(X1(i) + Y1(i)*tand(PA)) < tol )
        Ym = Y1(i);
        Xm = X1(i);
        i= length(X1)+10;
    end
end
flag = 0;
%Posizione angolare: non è nota a priori, ma sicuramente è compresa fra 0 e
%360/Z per periodicità dell'ingranamento 
teta = deg2rad(linspace(1, 360 / Z, 1000));
for k = 1 : length(teta)
  [XX2, YY2] = rotateZ(X2, Y2, teta(k));
  XX2 = XX2 + O2;
%   figure(2)
%   clf;
%   set(gca, 'color', [0 0.1 0.1])
%   set(gca, 'XColor', [0 0.5 0.5]);
%   set(gca, 'YColor', [0 0.5 0.5]);
%   axis equal;
%   grid on; hold on;
%   plot(X1, Y1, '-w')
%   plot(XX2, YY2, '-y')
%   plot(Xm, Ym, 'og')
%   plot([O2, -O2]*tand(PA), [-O2, O2], '-r')
%   xlim(.3*[O1, O2])
%   ylim(.5*[-O2, O2])
%   pause( 5 / length(teta))

    for i = 1 : length(X2)
       if abs(XX2(i) - Xm) < tol && abs(YY2(i) - Ym) < tol
           flag = 1; break
       end
    end

    if flag == 1, break, end
end
if flag == 0 ,  msgbox('CANNOT MESH','ERROR'), end

angle = teta(k);

end