function [Gx,Gy] = ColorGradients(I,t)

Ir = BiSmooth(I(:,:,1),4);
Ig = BiSmooth(I(:,:,2),4);
Ib = BiSmooth(I(:,:,3),4);
[rx,ry] = Gradient(Ir);
[gx,gy] = Gradient(Ig);
[bx,by] = Gradient(Ib);
cx = rx+gx+bx;
cy = ry+gy+by;

a = rx.^2 + gx.^2 + bx.^2;
b = rx.*ry + gx.*gy + bx.*by;
c = ry.^2 + gy.^2 + by.^2;

tr = a + c;

ind = find(tr > t);

Gx = zeros(size(rx));
Gy = zeros(size(rx));
for i=1:length(ind),
    S = [a(ind(i)) b(ind(i)); b(ind(i)) c(ind(i))];
    [V,D] = eig(S);
    if (D(2,2) > D(1,1))
        m = D(1,1);
        D(1,1) = D(2,2);
        D(2,2) = m;
        e = V(:,1);
        V(:,1) = V(:,2);
        V(:,2) = e;
    end;
    e = V(:,1)*sqrt(D(1,1));
    Gx(ind(i)) = e(1);
    Gy(ind(i)) = e(2);
    if (cx(ind(i))*Gx(ind(i))+cy(ind(i))*Gy(ind(i))<0)
        Gx(ind(i)) = -Gx(ind(i));
        Gy(ind(i)) = -Gy(ind(i));
    end;
end;

