I am trying to draw a tangent line to the bottom right of my circle and use the tangent line to find the angle between the horizontal line and the curve. Can someone please tell me what is wrong with my code.

% Contact Angle Detection

clc clear all close all

Irgb = imread('Trial.png');

% Processing of Real Life Image

Igray = rgb2gray(Irgb);

% Edge Detection Method

BW = edge(Igray,'canny');

% Finding Lines within the Image

[H,theta,rho] = hough(BW); imshow(imadjust(mat2gray(H)),[],... 'XData',theta,... 'YData',rho,... 'initialmagnification','fit'); colormap(gca,hot);

P = houghpeaks(H,5,'threshold',ceil(.7*max(H(:)))); lines = houghlines(BW,theta,rho,P,'Fillgap',5,'Minlength',7); imshow(BW), hold on; max_len = 0; for k = 1:length(lines) xy = [lines(k).point1; lines(k).point2]; plot(xy(:,1),xy(:,2),'linewidth',2,'color','green');

    plot(xy(1,1),xy(1,2),'x','linewidth',2,'color','yellow');
    plot(xy(2,1),xy(2,2),'x','linewidth',2,'color','red');
    len = norm(lines(k).point1 - lines(k).point2);
    if (len > max_len)
        max_len = len;
        xy_long = xy;
    end
end
plot(xy_long(:,1),xy_long(:,2),'Linewidth',2,'Color','Red');

x1 = xy_long(1); x2 = xy_long(2); x3 = xy(1); x4 = xy(2); y1 = xy_long(3); y2 = xy_long(4); y3 = xy(3); y4 = xy(4);

% Display Only the Curve Left Side rect = [300 323 50 10]; BW5 = imcrop(BW,rect); imshow(BW5); [y,x] = find(BW5);

% Polynomial Fitting of the Curve

p = fit(x,y,'poly3'); p1 = [3*p(1) 2*p(2) p(3)]; derivative = polyval(p1,unique(x)); avgDerivative = ((derivative(1)+derivative(2)+derivative(3))/3); CorrectedDeriv = -avgDerivative/10000;

imshow(BW) hold on line([x2 (x2-10)],[y2 (y2+(10*CorrectedDeriv))])

m2 = 0;

theta = atand(((CorrectedDeriv)-m2)/(1-(CorrectedDeriv)*m2));

%Display Only the Curve Right Side

rect = [463 323 50 10]; BW5 = imcrop(BW,rect); imshow(BW5); [c,d] = find(BW5);

% Polynomial Fitting of the Curve

p11 = fit(d,c,'poly3')

p12 = [3*p11(1) 2*p11(2) p11(3)] derivative1 = polyval(p12,unique(d)) avgDerivative1 = ((derivative1(1)+derivative1(2)+derivative1(3))/3) CorrectedDeriv1 = avgDerivative1/100

imshow(BW); hold on; line([x3 (x3+10)],[y3 (y3+(10*-CorrectedDeriv1))]);

m2 = 0;

theta1 = atand(((CorrectedDeriv1)-m2)/(1-(CorrectedDeriv1)*m2));

ThetaAvg = ((theta + theta1)/2)

[Image Analyst]

http://www.mathworks.com/matlabcentral/answers/13205#answer_18099