%
%    This file is meant to be used together with hf.m
%
%    It finds a loop between the generators pone and ptwo.
%
%    Ciprian Manolescu - May 18, 2005
%
%

Path = [];

% paths on the alpha curves going from pone to ptwo, left to right
% Path().alpha describes the vertices (as they appear in V)
% Path().alphaind describes the indices on the alpha curve
for i=1:g
	vect = [];
	wect = [];
	% ind describes the current point on the alpha curve 
	ind = V(G(pone, i), 2);
	endpt = V(G(ptwo, i), 2);
	j=1;
	while 1
		vect(j) = Av{i}(ind); 
		wect(j) = ind;
		if (endpt == ind) break, end
		if ind == (Alpha(i))
			ind = 1;
		else	
			ind=ind+1;
		end
		j = j+1;
	end
	Path(i).alpha = vect;
	Path(i).alphaind = wect;
end


% paths on the beta curves going from pone to ptwo
% Path().beta describes the vertices (as they appear in V)
% Path().betaind describes the indices on the beta curve

for i=1:g
        vect = [];
        wect = [];
        % ind describes the current point on the beta curve
        ind = V(Gbeta(pone, i), 4);
        endpt = V(Gbeta(ptwo, i), 4);
        j=1;
        while 1
                vect(j) = Bv{i}(ind);
                wect(j) = ind;
                if (endpt == ind) break, end
                if ind == (Sizeofbeta(i))
                        ind = 1;
                else
                        ind=ind+1;
                end
                j = j+1;
        end
        Path(i).beta = vect;
        Path(i).betaind = wect;
end


% the vectors LA and LG contain the intersection numbers between the 
% preliminary loop (described by Path) and the alpha and gamma curves, 
% respectively

LA = zeros (g,1);
LG = zeros (g,1);

% contribution from each piece of the beta curves
for i=1:g
	 for j=1:length(Path(i).beta)
		k = Path(i).betaind(j);
                if (Beta(i).curve(k)) > 0
                        ind = Beta(i).curve(k);
                        LA(ind) = LA(ind) + Beta(i).sign(k);
                else
                        ind = -Beta(i).curve(k);
                        LG(ind) = LG(ind) + Beta(i).sign(k);
                end
        end
end

% correction terms from the alpha curves
for i=1:g
	if any((Path(i).alphaind) == 1)
		LG (i) = LG(i) -1;
	end
	startpt = Path(i).alpha(1);
	whichbeta = V(startpt, 3);
	whichpoint = V(startpt, 4);
	startsign = Beta(whichbeta).sign(whichpoint);
	l = length(Path(i).alpha);
	endpt = Path(i).alpha(l);
	whichbeta = V(endpt, 3);
        whichpoint = V(endpt, 4);
        endsign = Beta(whichbeta).sign(whichpoint);
	if (endsign == startsign)
		LA (i) = LA(i) - endsign;	
	end
end

% decide if the generators are in the same spin^c structure

Corr = inv(BA') * LA;

if all ((abs(Corr - round(Corr))) < .0001)
	spinc = 1;
else 
	spinc = 0;
end